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1.0 Program ObjiAfives and Approach 

The primary goal of this research is to develop a solid-state high definition television 
(HDTV) imager chip operating at a frame rate of about 170 frames/sec at 2 Megapixels per 
frame. This imager offers an order of magnitude improvement in speed over CCD designs 
and will allow for monolithic imagers operating from the IR to the UV. 

The technical approach of the project focuses on the development of the three basic 
components of the imager and their integration. The imager chip can be divided into three 
distinct components: 1) image capture via an array of avalanche photodiodes (APDs), 2) 
charge collection, storage and overflow control via a charge transfer transistor device (CTD) 
and 3) charge readout via an array of acoustic charge transport (ACT) channels. The use of 
APDs allows for front end gain at low noise and low operating voltages while the ACT 
readout enables concomitant high speed and high charge transfer efficiency. Currently work 
is progressing towards the development of manufacturable designs for each of these 
component devices. 

In addition to the development of each of the three distinct components, work towards 
their integration is also progressing. The component designs are considered not only to meet 
individual specifications but to provide overall system level performance suitable for HDTV 
operation upon integration. The ultimate manufacturability and reliability of the chip 
constrains the design as well. 

The progress made during this period is described in detail in Sections 2-4. 
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2.0 Heterojunction Acoustic Charge Transport (HACT) Device 

In this section of the report we describe the progress made in the development of the 

ACT component of the imager chip. In general, this includes: 

1) The development and experimental investigation of new HACT architectures and 
the establishment of a Region of Operation curve, further, comprehensive. 

2) The completion of a comprehensive study of the piezoelectric thin film ZnO on 
Ga As-based materials for HACT applications. 

3) toitial studies into acoustoelectic amplifiers which are expected to drastically 
reduce the power budget for ACT devices. 

4) Advances in ACT device modeling which allow us to simulate device performance 
in any of the competing ACT architectures. 

2.1 New HACT Architectures and Established Region of Operation 

This section presents a physically-based analytic model for charge injection and 
transport in HACT devices. Previous models have utUized numerical solutions of the coupled 
Poisson and current continuity equations [2.1.I-2.1.3], This work, which was initially 
presented in [2.1.4], is the first physically-based analyUc description of the charge injection 
process for ACT devices. AnalyUc models for charge transport in buried channel ACT 
devices were investigated by Hoskins, e, al. [2.2.5]. This work extends this by refining some 
of the assumptions made by Hoskins, and by applying his approach to HACT devices. In the 
derivaUon of this model, bounds on different material parameters are developed that define the 
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range under which acoustic charge transport may occur. This model is verified 
experimentally. Because the conditions for charge transport must exist before injection can 
occur, this is treated first. Expressions for the equilibrium and nonequilibrium carrier 
concentration are then developed. From this, an expression for the instantaneous channel 
current is then derived. The actual shape of the injection aperture has been the subject of 
some controversy [2.1.6], [2.1.7], [2.1.1], [2.1.8]. It is found and verified experimentally that 
the shape is approximately Gaussian. The form of the steady-state channel current has also 
been debated [2.1.6], [2.1.1], [2.1.8]. It is shown that the steady-state channel current has an 
exponential dependence on gate voltage. 

The resulting models can be used to predict both the steady-state and transient 
characteristics of a HACT device of a particular geometry and material architecture. Being 
physically-based, generalizations about HACT device operation and trade-offs in device and 
material design can be made. Lastly, the model lends itself to a simple implementation on 
commercially available circuit analysis software. 

2.1.1 An Analytic Theory of Charge Transport in HACT Devices 

Hoskins and Hunsinger developed the first analytic theory of charge transport in ACT 
devices in [2.1.5]. This theory describes the transport of charge by the acoustically induced 
potential wave as a potential screening problem. In this section, this theory is extended and 
applied to heterostructure ACT devices. Different material constraints are identified, and their 
impact on the transport process is quantified. The model is also experimentally verified. 
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Figure 2.1.1 Potential screening in HACT devices. 
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2. 1.1.1 The HACT Potential Screening Problem. 

Acoustic charge transport is essentially a process whereby mobile electrons screen 
disturbances in the fixed charge concentration created by the motion of the SAW wave in the 
piezoelectric semiconductor [2.1.6]. Generally, this process is solved under the assumption of 
a stationary potential well, since the period of the SAW is much greater than the dielectric 
relaxation time of the medium [2.1.5]. The statement of the problem is then as follows: Given 
an acoustically- induced charge disturbance which is assumed stationary and superimposed on 
a background of fixed charge, determine the equilibrium free carrier concentration which screens 
this charge and the distance over which it is distributed. 

This problem is illustrated in Figure 2.1.1, which also shows the effects of the various 
layers in a typical HACT material structure. The primary purpose of the heterostructure is to 
confine mobile charge in the vertical direction [2.1.7]. A cap layer is used to control the effects 
of surface states on the conduction band in this vertical direction. Beneath this, is the charge 
control layer, which acts to form the upper potential barrier. Two alternatives have been 
proposed for this layer. Tanski, et al. proposed an Al 3 Ga 7 As layer which creates a potential 
discontinuity of approximately 0.25 eV in the conduction band [2.1.7]. Grudkowski, et al. later 
proposed a p-doped GaAs layer to create a barrier from the built-in potential of the depletion 
region [2.1.9]. Smith, et al. have also investigated this structure [2.1.10]. Ideally, the transport 
layer is the region where most of the mobile charge is confined. Beneath this, the barrier layer 
defines the lower region of charge confinement. A buffer layer is sometimes employed to 
minimize growth defects in the upper layers. 
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The potential to be screened is due to the superposition of the background donors and the 
polarization charge induced by the acoustic wave. The amplitude of the surface potential due 
to the acoustic wave, is related to the voltage, V„, applied to a transducer by 


where 11 is a parameter dependent on the metalization ratio (-0.8 for equal gap and width size). 
; is the electromechanical conversion efficiency (-0.7 typical). k= is the piezoelectric coupling 
constant (=0.0007 for GaAs) and N is the number of finger pairs in the IDT |2. 1 . 1 1 ] . The factor 
of two arises because the use of a reflector is assumed, which doubles the amplitude of the SAW. 


A common approximation assumes 


that the acoustic wave emanating from the IDT is an infinite 


plane wave in the transverse 


direction [2.1.12]. As a result of this approximation, the variation 


of the potential in the direction perpendicular to the transport direction is assumed to be zero. 
A similar, but somewhat cruder approximation is made with regard to the potential variation into 


the substrate. The amplitude of the potenual induced by a SAW on an electrically open surface 
is shown to fall off approximately exponentially with the acoustic wavelength into the substrate 
[2.1.1 1]. However, because the transport channel layer is usually a small fraction of an acoustic 
wavelength, the variation in potential in this direction can be neglected. Furthermore, the depth 
of the channel is generally small enough that the potenUal at the transport depth is essentially the 
same as at the surface, as given by equation (2.1.1). The SAW can therefore be expressed by 
an effective charge concentration, N,., related to the SAW potential by the one-dimensional 

Poisson’s equation 
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( 2 . 1 . 2 ) 



where is the acoustic wavelength, and is the peak polarization charge 

concentration, with q being the charge on an electron, and E, the permittivity of the material. 

The equilibrium electron concentration can be accurately obtained by self consistently 
solving the Poisson equation coupled to the Boltzmann distribution using numerical methods 
[2.1.13]. This has been done by a number of authors [2. 1.1-2. 1.3, 2.1.10]. However, insight into 
the physics of the device can also be gained by analyzing the asymptotic solutions of the 
nonlinear ordinary differential equation derived from directly substituting the Maxwell-Boltzmann 
distribution in Poisson’s equation [2.1.14]. The free electron concentration, n, can be expressed 
as a function of the conduction band potential \|/(x) = E^x)lq, where is the conduction band 
energy relative to the degeneracy level. 


n(x) = N^e 




( 2 . 1 . 3 ) 


where k is Boltzmann's constant, T is the absolute temperature, and N^- is the density of states in 
the conduction band. Assuming a transport layer doping, and neglecting any bandgap 
narrowing effects, equation (2.1.3) can be substituted in Poisson's equation, to yield 


dx^ 


P - Q 





-It 
- N,e " 

1 J 



( 2 . 1 . 4 ) 


where p is the charge density in the region of interest, equation (2. 1 .4) is a nonlinear ordinary 
differential equation that is not easily solved by analytic methods. To quantify the potential, and 


7 



hence the concentration throughout the transport channel, the domain over which equation (2. 1 .4) 
is solved is broken up into two distinct regions. Over regions where there is no gradient in 
potential (i.e. the "flat band" condition), equation (2.1.4) reduces to an algebraic expression, and 
the region is said to be quasi-neutral. This implies perfect screening of the SAW potential and 
donor ions by the free electrons. Hence, 


«(:«) = ^Dck 



( 2 . 1 . 5 ) 


This implies that the potential within the charge packet is 


kT 

i|t(a:) = — In 
<1 


N, 




V I 


0 ^ JC < .*0 


( 2 . 1 . 6 ) 


where defines the charge packet boundary. Hoskins, et al. arrived at the above result, but 
assumed that the potential outside the packet boundary was equal to the SAW potential [2.1.5]. 
This is true only for very low doped channels. Schmuckler later developed an empirical model 
for the problem of barrier lowering by the packet charge [2.1.15]. This work extends the 
previous work of others by developing a physically-based expression for the potential outside the 
charge packet by assuming that the depletion approximation holds whereby n is assumed 
negligible compared to the doping and polarization charge. In this case, the linear Poisson 

equation 
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(2.1.7) 




+ N^cas 



can be easily solved to arrive at an analytic expression for the potential. This can be 
accomplished by direct integration, yielding 


= + ^0 ^ ^ < V2 (2.1.8) 


where A and B are constants of integration. These cire evaluated by requiring continuity of the 
potential function and its derivative at Xg. In doing so. 



The only parameter left to determine to complete the approximate analytic solution to the 
potential is the length of the packet boundary Xq. Hoskins, et al. obtained a simple form of this 
screening length by assuming that the concentration drops to zero outside the region where the 
sum of the SAW potential and the ionized dopant potential goes negative [2. 1 .5]. Hoskins' model 
therefore assumes that charge can be screened, and hence transported in an undoped channel. 
In HACT devices, it has been experimentally determined that this is not the case. This will be 
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discussed in SeCion 2.1.3. A physical basis for this can be derived from Hutson and White's 
theory of acoustoelectric attenuators and amplifiers: Charge screening in HACT devices occurs 
only in regions where ihe spatial variaiion of the SAW potential is small on a Debye length scale 
within the channel [2,1.16], [2.1.17]. Thus, when the Debye length becomes an appreciable 
fraction of the acoustic wavelength, the polarization charge is not effectively screened by the free 
electrons in the channel. This forms a basis to determine v.. Debye length is a function of free 
carrier concentration and is defined by the length over which a small [«kT/„) perturbation from 
the flat bmid potenUal decays to 1/e of the initial value [2.1.14]. This results from the Taylor 
expansion of the exponential term in equation (2.1.3), and is not really appropriate for a large 
potential deviation. Considering this, a function, L„(x), termed the "quasi-Debye length," can be 
developed by substituting the packet concentration for the carrier concentration in the normal 

Debye length expression; 


in(4) 




e,kT 




(2h 


( 2 . 1 . 10 ) 


The screening length r, is defined by the condition L„(x„) - where is an empmcal 
parameter expected to be a small fraction of the acoustic wavelength X, In doing this, the 

screening length is given by 
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Figure 2.1.2 Predicted charge screehing length (packet length) as a function of SAW potential 
and channel doping. 


11 


( 2 . 1 . 11 ) 


K -1 


— cos 
2% 

q^LjjoN^ 


Figure 2.1.2 shows the predicted screening length (i.e. the length of the charge packet) versus 
SAW potential and channel doping for the material design presented in [2.1.10] . It is seen that 
there are local regions where the screening length is zero, indicating the charge transport will not 
occur. For this material design, the SAW potential threshold for charge transport to occur is on 
the order of 0.25 V in low doped channels. 

2. 1.1.2 Material Constraints and Operating Conditions for Efficient Transport. 

Eqs. (2.6, 2.8, 2.9 and 2. 1 1 ) form a complete analytic model for charge transport in HACT 
devices. The determination of L^o must be done empirically. The accuracy of this model is 
verified by comparing the predicted potential with that calculated by numerically solving equation 

(2.1.3) . This is shown in Figure 2.1.3. The Runge-Kutta method was used to solve equation 

(2.1.3) over one half of an acoustic wavelength [2.1.18]. It is seen that the analytic theory is 

close to that predicted by the numerical method. A value of L^o = V18 was found to give the 
closest fit for lower doped channels, while L,, = XJ5 gave a better fit for higher doped channels. 
It is also apparent that increasing the channel doping has two effects: increasing the well charge 
capacity, and decreasing the potential barrier between wells. Obviously, when the 

well is essentially full. An approximate expression for the full well charge capacity is easily 
derived by integrating equation (2.1.5) over the effective screening length given in equation 

( 2 . 1 . 11 ): 
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Potential vs. Distance 



(a) 



(b) 


Figure 2.1.3 (a) Calculated potential over one half wavelength for two different channel dopings, 
(b) Free electron concentration calculated from potential. 
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( 2 . 1 . 12 ) 



2h€ . 




\k\\ 


where the charge outside the screening length is neglected, and the concentration is assumed 
constant over the channel layer. The later assumption is only valid over the region where the 

channel is quasi-neutral. 

Hoskins, el al. pointed out that a trade-off must be made between the well capacity and 
the diffusion induced transfer inefficiency [2.1.5], The analytic model just developed provides 
a eonvenient means to predict the barrier height, and hence the CTE. Taking the barrier height 
as the difference between the potential minimum at jc=0, and maximum at x=XJ2. the well barrier 

height (|)„ is given by 




w 




( 2 . 1 . 13 ) 


A plot of the barrier height versus channel doping is shown in Figure 2.1.4 for several values of 
SAW potential. Hoskins used thermionic emission theory to predict the charge loss per clock 
cycle 12.1.5]. The results of his work indicate that barrier height of less than about 5kT may 
result in a significant loss of charge in typical device architectures. Figure 2.1.4 also points out 
several areas in which this theory improves upon models suggested m [2.1.5] and [2.1.6]. The 
maximum barrier height is somewhat lower than that predicted by these authors because of the 
inclusion of the effects of finite Debye length. For the same reason, the optimum channel doping 
is not equal to the polarization charge amplitude as they suggested. The optimum channel doping 
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Figure 2.1.4 Height of potential barrier between adjacent packets for several channel dopings. 
Polarization charge amplitude for a given SAW potential is inthcated by the squares. 
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would lie on a line passing through the peaks of the curves shown in Figure 2.1.4. This is seen 
to be somewhat less than the polarization charge amplitude. 

2.1.2 Experimental Technique for Transport Model Verirication 

The model developed in the last section was verified to some extent by comparing the 
predictions of the analytic expression to numerically derived results. However, to truly assess 
the validity of the model, the predicted values of several parameters must be compared to those 
obtained by measurements on actual devices. Basically, the transport model predicts the 
conditions under which acoustic charge transport will occur in HACT devices. It will be shown 
in this section that devices operating outside a range of SAW potential and equilibrium channel 
concentration are not measurably functional, whereas those that fall within this range show 
measurable ACT effects. Furthermore, equation (2.1.12) allows the prediction of charge capacity. 
This will also be verified experimentally. 

2. 1.2.1 Test Setups and Procedures. 

In order to verify that charge is being transported by the acoustic potential, two things 
must occur: 1) Channel current must be present only when the SAW potential is applied, and 

2) A transient signal imposed on the input must be delayed by the time it takes the acoustic wave 
to traverse the channel Figure 2.1.5 shows the setup used to measure both DC channel current 
and frequency response. While DC channel current is relatively easy to verify, the signal delay 
is sometimes not. Because of the low levels of the measured signal compared to spurious 
feedthrough signals, the expected frequency response due to the NDS array is badly distorted. 
However, it is possible to distinguish the acoustically delayed signal from the feedthrough signal 
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HOT DC TEST AND FREQUENCY RESPONSE SETUP 



Figure 2.1.5 Test setup used to measure DC channel current and frequency response. 
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by performing a Fourier nansform on .he measured fiequency response. This is sho»n in Figure 
2.1.6. It is seen that the feedthrough signal is of much greater amplitude dtan the acoustically 
delayed signal. However, if the output taps are spaced sufFtciently far from rite input gate, the 
feedthrough signal decays to the point where the acoustically delayed signal is measurable. 

Hence the ACT operation is verified. 

2. 1.2.2 Measurements and Results. 

Figure 2.1.7 shows a summary of the performance of HACT devices fabricated a. Georgia 
Tech and United Technologies versus the equilibrium channel concentration and SAW potential. 
This plot shows three different consriains on the channel concenttation and SAW potential. Two 
of these, the Polarization Charge Limit and the Transport Limit are discussed in the section, while 
the injection Barrier Limit is discussed in Section 2.1.3. A lower limit on charge nansport is 
defined by rite screening length, as given in equation (2.1.11), reducing to zero. Since this is 
dependent on frequency, the three different curves represent the diree different operating 
ftequencies of the different devices measured. Filled markers indicate device wafer lots in which 
charge transport was verified, and empty markers designate nonfunctional devices. 1. is seen that 
this limit is only relevant at low SAW potentials « 8kW«), and only for low channel dopings. 
It is seen that the condition is automatically met if the conditions are met for Schottky injection, 
which is discussed in the next section. Thus it was not possible to verify the lower boundary of 

HACT operation using the injection stnictures of the test devices. 

VerificaUon of the full-well charge capacity given by equation (2.1.12) was done by 
n^asuting the maximum channel current that the SAW can cany . Figure 2.1.8 shows the channel 
current as a function of gate voltage for United Technologies HACT B-3. This will be discussed 
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Figure 2.1.6 Output signal of device D02 of HACT lot 1854. Time domain signal is calculated 
by a Fourier t rans form of measured frequency response. ^Signal delay of 600 ns confirms ACT 
operation 
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l..r. 2 1 7 Functionality of HACT devices veisiis SAW potential and ch^el concentrations 
iUed miters indicate that charge transport was verified, empty markers denote non-func 
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in gn^atnr dclail in the next section. It is seen that the channel current abruptly stops increasing 
as a function of gate voltage, indicating that the potential wells are at their peak capacity for the 
operating SAW potential. Since the average channel current is simply the product of the packet 
charge and the SAW frequency, the right axis is plotted in units of electrons per packet. It is 
seen that the maximum capacity is close to that ptedicted by equation (2.1.12). 


2.1.3 Development of the Injection Model 
As illustrated in Figure 2. 1 .9, the input structure superficially resembles a HEMT [2.1.19]. 
However, differences between the two devices are numerous. Under normal operation all charge 
is moved by the SAW potential travelling at the acoustic velocity, which is two orders of 
magnitude slower than the saturation velocity in a HEMT. The epitaxial layers also differ. In 
a HEMT, the AlGaAs charge control layer in heavily doped n-type. As described by the 
Anderson model, the differing electron affinities of the narrow bandgap GaAs and the wider gap 
AlGaAs create a conduction band discontinuity at the heterojunction interface [2.1.20]. In 
equilibrium the doped AlGaAs layer provides electrons which collect in the potential well formed 
by this discontinuity. This creates a highly concentrated 2-DEG channel region. In contrast, the 
charge control layer in a HACT acts only to form a space-charge region which broadly confines 
the carriers in the vertical direction the transport channel. Also, as discussed in the last section, 
the equilihrium channel concentration is kept quite low (-10“ cm > ) so as not to exceed the 
polarization charge limit. Moreover, the cap layer of a HEMT structure is heavily doped n-type 
to aid in the formation of channel ohmic contacts. This is not possible in the HACT stmcture, 
because the shorting of the acoustic potential at the surface will cause a significant decrease in 
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Figure 2.1.9 Schematic of the HACT injection structnre and associated equilibrium band 
diagram. 
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.he SAW potential at the transport depth 12.1.11). Thus it turns out that most MESFET and 
HEMT equivalent circuit models are no. really applicable to HACT devices. However, the 
subthmshold region of MESFET and HEMT operation, no. generally considered in the 
development of equivalent circuit models, does resemble die charge injection into the HACT 
channel. When the channel region is completely pinched-oft. die MESFET or HEMT is said to 
be in the mbthrer/toW region (2.1.14]. [2.1.211.(2.1.22). In a conventional MESFET or flEMT 
operating in this regime, charge present in the n* ohmic region diffuses over the potential barrier 
set up by the gate voltage. It might also be noted that injection in HACT devices is similar to 

charge injection in heterostructure insulated gate FETs (HIGFET's) (2.1.23). 

The theory of HACT charge injection is presented in the following subsections. The 
development diis dieory parallels the theory of subdueshold injection in more conventional 
semiconductor devices. The inclusion of the SAW potential does complicate the problem 
somewhat, so the derivation is broken up into several steps, the fimt of which is to determine 
the equilibrium mobile electron concentration duoughout die injection region. Although 
numerical methods are available to do diis accurately, some assumptions can be made that allow 
the formulation of an approximate analytical solution. The next step is to evaluate the 
nonequilibrium concentrations in the injection region. Again, simplifying assumptions are made 
to make the problem tractable using analytic medtods. Once the mobile carrier concentrations 
are known, the final step is to incorporate these into the drift-diffusion equation in order to arrive 
a. an analytical expression for the instantaneous channel current. An additional step is performed 
to obtain a closed-form expression for the sready-sute channel cureenl as a function of gate 

potential. 
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2. 1.3.1 Determination of Equilibrium Carrier Concentration. 

Referring to Figure 2.1.6, the injection region is divided into four parts. In the first 
region, the ohmic contact (xkLq^) contains a high concentration of mobile electrons due to the 
diffusion of dopants during the alloy process. The determination of the actual free carrier 
concentration in this region is quite difficult for several reasons. The first of which is the 
degenerate nature of the doping in this region Nds(x), which precludes the use of Boltzmann 
statistics. Secondly, widely differing reported experimental results imply that the alloy process 
itself is not well understood [2.1.24]. It is often assumed that the concentration is on the order 
of 10” cm'^, and that the ohmic drop during the current injection is negligible [2.1.25]. Within 
the ohmic region, the conduction band is considered flat, and quasi-neutrality is assumed. Hence, 
H = since the doping concentration is many orders of magnitude higher than the polarization 
charge. The potential is related to the carrier concentration by 


nix) = 



( 2 . 1 . 14 ) 


where is the Fermi integral of order 1/2 [2.1.14]. 

The second region, between the source and gate contacts ( < Jt < 0), will be referred 

to as the "injection region" throughout this section, as it is the concentration in this region that 
most greatly affects the charge injection operation. Here the concentration can be found from 
the solution of Poisson's equation. Since the potential in this region is reacting to a large change 
in doping, neither the quasi-neutrality nor depletion approximations can be assumed [2. 1 . 14]. The 
equation to be solved is then 
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(2.1.15) 


d:? 


Noix) - N^e 


where = N.. for . < -L,, and = N,.* for -L,, < . < 0. What is often done for small 
potential perturbations {A^v«kT/q) is to expand the exponential dependence of the carrier 
concentration on the potential in a Taylor series, resulting in a linear ordinary differential 
equation, and the well-known definition of Debye length [2.1.19]. This result can be applied to 
the more general problem by recognizing that a larger change will result in a similar exponential 
change in the potential. The form of this dependence can be deduced by empirically altering the 
solution to the linear problem as follows: 


where 
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(2.1.17) 


Requiring continuity of the potential and its derivative results in expressions for the parameters 
a, and ot 2 '. 
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(2.1.18) 
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where y, and y^ are determined empirically. One might expect that the potential will drop off 
somewhat less than the Debye length scale in the source ohmic region, and somewhat greater in 
the injection region. This implies that yi<l, while y 2 >l. Figure 2.1.10a shows the potential 
calculated from equation (2.1.16), and that calculated from PC- ID, a numerically-based Poisson 
solver [2.1.26], The x-axis is scaled by It is seen that for y,=0.05 and y 2 = 8 , good 

agreement is obtained over the injection region. The larger discrepancies around Jt=-Lgj result 
from degenerate effects not included in equation (2.1.16). The concentration in the injection 
region is obtained from substituting equation (2.1.16) in the Boltzmann distribution. 


n(jc) = fV^exp 
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(2,1.19) 
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Figure 2.1.10b shows the concentration as calculated from equation (2.1.20) compared to that 
calculated using PC-ID. The agreement is good over most of the injection region. 

The next region to discuss is that beneath the gate (OcrcLg). Under normal operation, 
this area is depleted by the built-in voltage of the Schottky contact. This forms the basis to 
derive the threshold voltage V,* by integrating Poisson’s equation. Variation of the potential in 
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(a) 



Distance pebye Lengths) 

(b) 


Fi*.;. 2.1.10 (,)Inj«Uon region 0» Electron 

concentration calculated analytically and from PC-ID. 




the jt-direction can be neglected if the gate is longer than a few Debye lengths, and the 
integration done in one dimension. In the absence of a gate bias, the total potential in this region 
is the superposition of the Schottky barrier, depletion, and conduction bandgap potentials, 
provided they exist. The threshold voltage then is given by the potential at the bottom of the 
depletion region, which is assumed to extend to the bottom of the channel layer. In the case of 
the p-HACT structure proposed by Smith, et al. [2.1.10], a pn junction is used to form the top 
barrier, and the gate is recessed to the charge control layer. The threshold voltage for this is 
given by 



where is the Schottky barrier height, is the thickness of the charge control layer, and 
its donor concentration, and is the thickness of the transp)ort channel layer with doping 
In the case where an AlGaAs layer is used to form the barrier at the top of the channel, the 
conduction band discontinuity AE^ must be included in the calculation of threshold voltage. The 
gate is not usually recessed in this material structure. Thus V,^ is given by 

^ ~^^c) * '^{^<P^Dcp ^cc^Dcc ^ck^Dch) ( 2 . 1 . 22 ) 

where d^^ is the cap layer thickness, and is its doping. 

The last region to consider is the channel region, x>Lq. The equilibrium concentration 
is obviously at a location several Debye lengths from the end of the gate. Being many 
orders of magnitude smaller than this creates a barrier in equilibrium <|)^, which must be 
overcome by the SAW potential. The height of this injection barrier is given by 
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( 2 . 1 . 23 ) 



where is the inverse of the Fermi integral function 




( 2 . 1 . 24 ) 


This is not an easily invertible function; however the tabulated values indicate that the conduction 
energy lies about 6kT below the degenerate level for doping on the order of 10 cm [2. 1.19]. 
As a minimum, the SAW potential amplitude (j),o must exceed the injection barrier (-0.35 V for 
cm \ or no net current will flow into the channel. The injection barrier is plotted as 

an operational boundary in Figure 2.1.7 in Section 2.1.2. 

2. 1.3. 2 Nonequilibrium Carrier Concentration. 

When a SAW potential (|>,(a, 0 is propagating through the HACT mesa, the concentration 

in the source and injection regions are only perturbed a small amount since the polarization 
charge concentration is many orders of magnitude smaller. Hence, the effects of the SAW 
potential are neglected here. As seen from the quasi^equiUbrium ACT model developed m 
section 2.1.1, the concentration in the channel region is greatly affected. A similar situation 
occurs in the gate region. When a gate bias voltage V<„ is connected between the source and gate 
contacts, a depletion depth changes beneath the gate. The total potential beneath the gate is now 
given by V„ - f .), assuming negligible current flow into the gate and negligible 
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screening from any mobile carriers that might be present. The free electron concentration at the 
input of the gate region, n(x=0,t) is then given by 
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exp 
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(2.1.25) 


assuming that L^s is small enough that the concentration is much greater than the channel doping 
and the acoustic polarization charge It is seen that the gate voltage controls the amount 
of charge available for transport by the SAW potential, if the conditions specified in section 
2.1.1 are met. 

2. 1.3.3 Drift-Diffusion Injection Theory. 

Having established the concentration in the injection region, the next step toward 
formulating the injection model is to analyze the motion of the electrons in the presence of the 
SAW potential. The drift-diffusion model under steady state assumptions is adequate to describe 
the motion of carriers in the injection region because of the slow variatic i of the piotential in the 
transport direction [2.1.5]. When a SAW potential with an amplitude in excess of the injection 
barrier is propagating, channel current is expected to flow from source to drain due to the carrier 
bunching effects discussed in the last section. The application of a negative gate bias will limit 
this current by deepening the gate depletion region, and thus heightening the barrier over which 
the carriers must flow. A common technique used to evaluate carrier transport over a barrier is 
to integrate the current over the barrier region. Since the current is constant, it can then be 
evaluated from the boundary conditions at the depletion region edges [2.1.14]. Assuming that 
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the carrier concentration in the GaAs transport channel is uniform in the j-direction, the 
instantaneous channel current I, through the gate region is 




I 

-Wd,,q\n(x,t)y^,— 



( 2 . 1 . 26 ) 


where W is the channel width, is the electron mobility, and D„ is the diffusion constant for 
electrons. In keeping with the development of the MESFET subthreshold model of Liang, er ai, 
the drift current term in equation (2. 1 .26) is neglected because of the low concentration of carriers 
in the depletion region, and because the potential beneath the gate is relatively flat [2.1.22]. 
Assuming that recombination in the depletion region can be neglected, the current continuity 
equation forces the current to be constant throughout this region. Hence, integration of equation 
(2.1.26) results in a simple expression for the injected channel current. 


An expression for the concentration at x=0 is given in equation (2.1.25). Initially, the well into 
which charge is injected into is empty, and n(i.„,0) « n(0,0). Hence, the channel cutrent can 


be approximated by 
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In reality, the well is filling up as charge is injected, hence decreasing the rate of charge 
diffusion. /, will decrease substantially from this value as the charge packet concentration 
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increases to a significant fraction of n{0,t). Hence, equation (2.1.28) is valid in the low injection 
limit. The rate equation can be solved to obtain a more exact form of the injected charge. 
However, this rate equation is nonlinear because concentration in the well changes as a function 
of the charge load. In order to keep a closed-form expression, a simple empirical modification 
to equation (2.1.28) can be made by introducing an ideality factor Ti: 


Wd,.qD^ iVnciVn,* 
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•exp 




(2.1.29) 


Strictly speaking, the ideality factor T| is a function of gate voltage through its relation to 
the well charge. However, proper selection of this parameter results in agreement over wide 
ranges of operation. It should be pointed out that some ideality factor would have to be 
introduced into equation (2. 1 .28) regardless of the well filling effect, to account for recombination 
in the depletion region. It has been found for depletion regions of lower doped pn junctions that 
T|=2 results in a good fit to the measured I-V characteristics because of the bulk recombination 
effect [44], So it is expected that ti=2 would be the asymptotic limit for low injection, and that 
T) would increase as the average injection current increased. 

2. 1.3.4 Gaussian Approximation for the Instantaneous Channel Current. 

The expression given in equation (2.1.29) for the instantaneous channel current is closed- 
form, but is not easily integrable because of the cosine term in the exponential. For the purposes 
of calculating the steady-state channel current, the cosine is expanded out to two terms in a 
Taylor series. A sampling aperture function S(t) can then be defined which is essentially 
Gaussian in form. 
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Normalized Time (Ts) 


Fi,u« 2.J.I1 Calculated sampliug. aperture fimcUon using exact fonnulation and CSaussian 
approximation for a 1 mm wide channel, r\ = 5.5 and <>, - 0.26 V. 
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( 2 . 1 . 30 ) 


where 



J _ ^GS 
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S(0 . e o = J- -!l*I . 

»/.N »♦« 

Plots of 7/r; calculated from Eqs. (2.24) and (2.25) are shown in Figure 2.1.1 1. It is seen that 
there is little error introduced by the Gaussian approximation. Defining the sampling width to 
be a as given above would imply that 67% of the charge injected into the well occurs during this 
time. The approximate width of the sampling aperture has been the subject of much debate. It 
has been claimed that the aperture is a very small fraction of a SAW period [2.1.6], [2.1.7], 
[2.1.27], From equation (2.1.31), typical values for HACT device designs reported in [2.1.7] 
result in apertures approximately 30% of the SAW period. This is experimentally verified in the 
next section. 

2. 1.3.5 Average Channel Current. 

The average channel current )) can be obtained by integrating over one SAW 

period . This can be done by numerically integrating equation (2.1.29), or analytically by 
expanding the limits of integration to infinity on the Gaussian definite integral formed from 
equation (2.1.30). In doing the later, an approximate analytical expression for the average 
channel current is developed; 
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equation (2.1.32) predicts an exponential dependence of the average channel current on the gate 
voltage, which is the same as observed in MESFETs and HEMTs in the subthreshold region. 
This is verified experimentally in the next section. 


2.1.4 Experimental Technique for Injection Model Verification. 

Verification of the analytical injection model presented in this chapter was made by 
measurement of HACT devices fabricated at Georgia Tech and United Technologies, Inc. The 
two salient features predicted by the injection model are the Gaussian shape of the sampling 
aperture and the exponential dependence of the average channel current on the gate voltage. 

2. 1.4.1 Injection Aperture Test Setups and Procedures. 

The Gaussian shape of the aperture function has been verified by measuring the 
conversion loss of the input sampler with input signal frequencies at harmonics of the Nyquist 
frequency. The mathematical basis of this measurement, and how measurement of the conversion 
loss relates to the aperture shape, can be understood by referring to Figure 2.1.12. This shows 
the functional block diagram of the HACT devices on which these measurements were performed. 
The output current is the product of the sampling function and the input voltage signal. Since 
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Figure 2.1.12 System level representation of the HACT device for aperture function 
measurement. 
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the sampling function is periodic, it can be expressed as a Fourier series 


Sit) = E 

m*-“ 

where = 2 ^ 1 / 5 , and a„ is the mth Fourier coefficient. Thus, the determination of 

these coefficients leads to the shape of the sampling aperture. Assuming that 5(0 is real and even, 
then a„>0 for all m. If a sinusoidal signal, v*(0, of a frequency some integer multiple of the 

Nyquist frequency 


v^(0 = cos 



( 2 . 1 . 34 ) 


is applied to the gate, then the output channel current, bandlimitted by the bandpass filter function 
of the NDS array, can be expressed as 

f.«) = 

where 5 is the Kronecker delta function. It is seen that the amplitude of iy is twice the value of 
the Fourier coefficient a„.y. Thus, measurement of the power at the kth multiple of the Nyquist 
frequency will give the Fourier coefficient, and hence the shape of the injection aperture. 

2. 1.4 .2 Injection Aperture Measurements and Results. 

The injection aperture was measured on two HACT devices with identical input structures 
using the setup shown in Figure 2.1.13. The time-domain waveform was reconstructed from 
measurement of the first three Fourier coefficients. The resulting aperture function as calculated 
using equation (2.1.33) is shown in Figure 2.1.14, along with the Fourier reconstructed sampling 

function. 
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Figure 2.1.13 Test setup used to measure injection aperture shape. Input signals are multiples 
of the Nyquist frequency. Output is input frequency translated by a harmonic of the aperture 
function. 
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2.J.14 Measured versus predicted aperture function for wo 1 mm wide HACT devices 



2. 1.4.3 Average Channel Current Test Setups and Procedures. 

The relationship between steady-state channel current and gate voltage was also verified. 
This was done using the test setup discussed in section 2.1.2. The average channel current is 
calculated from the approximate Gaussian integral given in equation (2. 1 .32), and by numerically 
integrating equation (2.1.30). The difference between these is too small to discern, hence both 
results appear as the "predicted" values in Figure 2.1.15. The measured curve of a 1 mm 

wide HACT device designed and fabricated by United Technologies is also shown in Figure 
2.1.15. The ideality factor T| and was adjusted so that the calculated (I) was in close agreement 
with the measured average channel current. It was found that T|=3.3 gave good results over a 
wide range of channel current values. 

2. 1.4.4 Average Channel Current Measurements and Results. 

The nature of the average channel current as a function of gate-source voltage has been 
the subject of some debate. An exponential dependence was shown by the numerical model 
developed by Bogus for buried channel ACT devices [2.1.1]. However, Merrit observed a 
square-law dependence in HACT devices, similar to what is seen in a MESFET above threshold 
[2.1.8]. It was observed in this research that the measured channel current does follow an 
exponential dependence as predicted by the injection model for low injection levels. However, 
the variation of the ideality factor with gate voltage causes a less than exponential dependence 
in the (O-Vcs characteristic. It is possible that this variation would result in an empirical fit to 
a square-law characteristic. However, the physical basis for this is due to completely different 
mechanisms than are present in a MESFET operating in the saturated region. 
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2.1.5 Summary of Results and Conclusions 


This chapter has presented the first analytic theory of charge transport and injection in 
HACT devices developed from the physics of the device operation. The transport model includes 
expressions for charge packet concentration, charge capacity, and barrier height between adjacent 
packets, in terms of the physical device parameters and operating conditions. It was seen that 
charge transport is greatly affected by channel doping, and there exist maximum and minimum 
conditions for both the doping and SAW potential in order for charge transport to occur. The 
theory was verified against measured data for a number of HACT devices fabricated at Georgia 
Tech and elsewhere. The injection theory predicts the instantaneous and average channel current 
as a function of device geometry and operating conditions. This was the first analytic description 
of the injection aperture shape, which was found to be approximately Gaussian. The physical 
mechanism of subthreshold current injection was also used to explain the average channel current 
dependance on gate voltage being close to exponential. Both the instantaneous and time-averaged 
channel current functions were verified experimentally. 

The work presented in this chapter provides a framework within which one can quantify 
tradeoffs in the design of HACT material architecture and device geometry. The model can also 
be used to predict the limits of the technology for specific applications. 
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2.2 Piezoelectric Thin Films on GaAs-Based materials 

A potential application for piezoelectric films on GaAs substrates is the monolithic 
integration of surface acoustic wave (SAW) devices with GaAs electronics. Knowledge of the 
SAW properties of the layered structure is critical for the optimum and accurate design of 
such devices. The acoustic properties of ZnO films sputtered on (001) -cut 
<1 10>-propagating GaAs substrates are investigated in this article, including SAW velocity, 
effective piezoelectric coupling constant, propagation loss, diffraction, velocity surface, and 
reflectivity of shorted and open metallic gratings. The measurements of these essential SAW 
properties for the frequency range between 180 and 360 MHz have been performed using a 
knife-edge laser probe for film thicknesses over the range of 1.6 to 4 pm and with films of 
different grain sizes. The high quality of dc triode sputtered films was observed as evidenced 
by high and low attenuation. The measurements of the velocity surface, which directly 
affects the SAW diffraction, on the bare and metalized ZnO on SiO^ or SijN 4 on {001 }-cut 
GaAs samples are reported using two different techniques: 1) knife-edge laser probe, 2) 
line-focus-beam scanning acoustic microscope. It was found that ne^ the <11 0> propagation 
direction, the focussing SAW property of the bare GaAs changes into a non-focussing one for 
the layered structure, but a reversed phenomenon exists near the <100> direction. 
Furthermore, to some extent the diffraction of the substrate can be controlled with the film 
thickness. The reflectivity of shorted and open gratings are also analyzed and measured. Zero 
reflectivity is observed for a shorted grating. There is good agreement between the measured 
data and theoretical values. 
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2.2.1 Background: 


Over the past two decades, ZnO has been the most frequently considered piezoelectric 
film for use for increasing the piezoelectric coupling in cases of a weakly- (e.g. GaAs) or 
non-piezoelectric substrate [2.2.1], The use of a ZnO layer with a GaAs substrate will 
enable one to monolithically integrate surface acoustic wave (SAW) devices with GaAs 
electronics. Furthermore, such a structure will have low dispersion due to the similarity of the 
SAW velocities of both materials [2.2.2]. In this paper, the feasibility of ZnO films deposited 
on GaAs substrates for SAW device applications is investigated along with a detailed 
discussion about the problems which occurred during the characterization of the ZnO/GaAs 
structures. 

Knowledge of SAW properties of the filmed structure is crucial for the accurate design 
of SAW devices. Herein, the experimentally determined fundamental SAW properties are 
reported including velocity, effective piezoelectric coupling constant K^, attenuation for the 
1.6 to 4 pm thicknesses of c-oriented ZnO film over semi-insulating [OOlJ-cut 
<1 10>-propagating GaAs substrates in the frequency range of 180-360 MHz. 

Another critical acoustic property for the design of SAW devices is diffraction of the 
surface acoustic wave. Diffraction in anisotropic media can generally be described by the 
slowness surface or the velocity surface. The velocity surfaces on free and metalized surfaces 
were measured using two independent techniques: 1) knife-edge laser probe, 2) 
line-focus-beam scanning acoustic microscope (LFBSAM). The results from both methods 
are compared with the theoretical results. 
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Metallic gratings are basic elements required for the construction of SAW devices. 
Analyzing the reflectivity and the velocity change due to metallic gratings has been studied 
since the invention of SAW devices; however, no research on gratings on multilayered 
structures has been reported. Dana and Hunsinger have analyzed the grating properties on a 
single crystal by separating into the piezoelectric shorting and the mechanical scatterrng due 
to the mass loading and the stress [2.2.3]. In this paper, an eztension of their technique to a 
multi-layered structure is presented. Experimental results for a shorted and an open grattng 
are compared with theoretical calculations. 

A passivation layer such as SiO; [2.2.4], SijN 4 [2.2.5], or SiON [2.2.6] between ZnO 
and GaAs might be required in order to: 1) passivate the structure and enhance the yield, and 
2) ptevent unwanted doping of GaAs by Zn during the sputtering processes. The application 
of such a passivation layer was obtained with a plasma-enhanced chemical vapor deposition 
(PECVD) for SiOj and CVD for SijN,. The acoustic parameter differences with and without 
the passivation layer are discussed. Theoretical calculations of SAW velocity and K’ as a 
function of ftlm thickness were performed using the Laguerre polynomial technique [2.2.7]. 

2.2,2 ZnO film growth 

C-axis oriented polycrystalline ZnO fllms were grown at Motorola using both rf 
magnetron and dc triode sputtering method on [001]-cut GaAs wafers which were grown by 
the liquid-encapsulated Czochralski (EEC) technique. The GaAs wafers, which were 
chemically and mechanically poUshed on both sides, were semi-insulating with resistivity d = 
lO’n-em, and the surface normal direction was [100] ± .1“. The quality of the ftlm suongly 
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depends on the fabrication conditions. The typical deposition parameters used for the ZnO 
film are listed in Table 2.2.1. The film thicknesses of 1.6, 2.8, emd 4.0 |im have been 
chosen to cover the range of 0. 1-0.5 of the acoustic wavelength, X. As a passivation layer, a 
lOOOA PECVD S 1 O 2 layer was deposited for the rf magnetron sputtered ZnO films and a 
2 OOOA CVD Si 3 N 4 layer was deposited for the dc triode sputtered films. The choice of a 
passivation layer was made solely on the equipment availability. 

The as-grown ZnO films were transparent with a very smooth surface finish. 

The grain size of the ZnO films was 0.2 to 0.5 pm for the rf magnetron sputtered films, and 
no grain boundaries are visible under scanning electron microscope (SEM) examination for 
the dc triode sputtered films indicating the grain size was beyond the resolving power of the 
SEM. From X-ray diffraction and atomic force microscope (AFM) measurements, the grain 
size of the dc triode sputtered films is 30 to 50 nm. The problem with these particular dc 
triode sputtered films, however, is that a compressive stress is induced in the substrate at the 
time the film is grown. The induced stress is greater with the thicker films such that 2" wafers 
are warped, making it difficult to process IDTs. Thus, only the ZnO film thickness of 1.6 pm 
was investigated for the dc triode sputtered films. It has been previously rejjorted that a 
decrease in the internal compressive stress results in a decrease in [2.2.8]. 

The orientation of the ZnO film was characterized using the X-ray diffraction 
technique. The diffraction patterns showed that all samples had strong c-axis texture; i.e., 
there were no ZnO peaks visible except those from the basal planes. 


■ 47 



2.2.3 Characterization of acoustic properties 
A number of devices were fabricated using a standard lift-off process on 1.6, 2.8, and 
4 pm thicknesses of the ZnO films. Employing one set of four IDTs with different center 
frequencies enables the investigation of frequency characteristics over the range of 180-360 
MHz. The configuration of the IDTs is listed in Table 2.2.2. The values of the SAW 
wavelength, X, in the table denote nominal values. Each IDT is unanodized with metallization 
ratio (finger width to grating period) of 0.6. The apertures of the IDTs were chosen to 
faciUtate the measurement of the slowness surface on hare GaAs substrates 12.2.9). 
Considering tbe non-acoustic resistivity of the IDTs (2.2.10), the metal thickness was chosen 
to be 1000 A of Al/4% Cu alloy. For ease of fabrication, the smallest dimension for the 
devices was no less than 2 pm; therefore, IDT No.l utilized spht finger electrodes while the 
others utilized the single finger configuration. Acoustic absorber was applied at both ends of 
the device in order to prevent SAW reflections from the edges. 

2.2.3. 1 Vsaw 

The value of X is inferred from the periodicity of the IDT electrodes, and the SAW 
velocity. in the IDT region may be obtained to first order by the relation, = ij'k, 
ivhere f, represents the center frequency. The value of the nominal X within the free 
ZnO/GaAs surface region is different and is measured directly by using a technique such as 
the laser probe system. The measured values of f, are shown in the Table 2.2.2, and the 
variation of f, is due to the fact that is dispersive depending upon the ZnO film 
thickness as well as the passivation layer. 
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A knife-edge laser probe system was used to measure SAW properties. The use of the 
laser probe system and the scheme for the longitudinal and transverse scans has been 
described in [2.2.1 1] and [2.2.12] and will not be discussed here. In this study, the signal to 
noise ratio of the laser probe was typically 40-50 dB, and a 3 p.m diameter laser spot was 
used. From several devices on each subs trate, the best device was chosen for the experiments. 
The SAW velocity measured with the knife-edge laser probe is shown in Fig. 2.2. 1 along with 
the theoretical values. The velocity dispersion is plotted as a function of 
wavelength-normalized film thickness. The theoretical values were calculated with the value 
of X set to 10 Jim. The material constants in the calculation were those given by Slobodnik 
[2.2.13] and Hickernell [2.2.5]. 

The maximum difference between the SAW velocities of the free surface ' measured 
with the laser probe) and in the IDT area, from the center frequency measurements, was 
0.3%. The maximum velocity deviations caused by the passivation layer were 0.6 and 1.6% 
for SiOj and Si 3 N 4 , respectively, and these figures agree well with the values from the center 
frequency measurements. The accuracy of the velocity measurement is approximately 0.3%. 
The detailed discussion about the accuracy can be found in Bright, Kim, and Hunt [2.2.14] 
and will not be reproduced here. 

In using the laser probe, difficulty in making the free surface measurements was 
experienced due to the multiple reflections of the laser beam between the film-substrate 
interface and the air-film. For three film structures (2.8 and 4.0 jxm rf magnetron sputtered 
ZnO and 1.6 ^im dc triode sputtered ZnO/0.2p,m Si 3 N 4 ), no appreciable energy from the laser 
beam was reflected from the free surface areas which is similar to what happens with 
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Figure 2.2.1 SAW velocity measured by a knife-edge laser probe on the metalized surface: 

(a) RF magnetron sputtered film, 

(b) DC triode sputtered film. 
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anti -reflection coatings in optics. This made probing impossible. For this reason, the velocity 
shown in Fig. 2.2.1 is the measured data on the metalized surface as well as the theoretical 
value on short-circuit surface. More detailed analysis about this multiple reflection can be 
found in [2.2.15]. 

2.2.3.2 Piezoelectric Coupling (K^) 

The effective piezoelectric coupling constant is generally obtained by a 
measurement of the relative shift in velocity between open- and massless short-circuit surface 
as follows, 

(2.2.1) 

where and are the open- and short-circuit velocity, respectively. For the measurement of 
V^, a metal pad was fabricated on one side of the IDTs to short out the electric field on the 
surface while the other side of the IDT was left as a free surface. In order to reduce the mass 
loading of the metal pad, it was fabricated with 250 A thickness of Al/4% Cu alloy. In spite 
of the finite thickness of the pad, its mass loading would be negligible because its thickness is 
only 0.0016 to 0.0031 of X (16 to 8 fxm) and it has a very slight effect compared to the 
piezoelectric effect [2.2.11]. This assumption was confirmed by a theoretical calculation using 
the Laguerre polynomial technique [2.2.7]. The velocity change due to the mass loading was 
calculated to be at most 0.015% which is much less than the expected value of for the 
ZnO/GaAs structure. 

The value of calculated by using Eq. (2.2.1) is shown in Fig. 2.2.2 as a function of 
normalized film thickness. The data points represent average values obtained from several 
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Figure 2.2.2 Piezoelectric coupling constant 
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measurements. It is seen that the value of for the dc triode sputtered film is about 1.7 
times larger than that of the rf magnetron sputtered films as expected from the grain size of 
the film. The measured value is 40-80 % of the theoretical one, which is comparable to 
Hickeraell's observations on the other substrates [2.2.1]. This value is still approximately 
5.7-10.8 times larger than the typical value for bare GaAs (0.07 %). 

2.2.3. 3 Propagation Loss 

Acoustic propagation loss (attenuation) is also a crucial parameter in device 
applications of piezoelectric films. The attenuation of the SAW can be measured by 
comparing the energy contained in the wave of two different transverse scans, separated by 
some distance d The energy is proportional to the value obtained by integrating the squared 
magnitude of the beam profile lUP along the transverse direction. The scan length should be 
made long enough to include the entire energy since the beam spreads out due to diffraction. 
Thus, the attenuation a is given by 


a=101ogj 


[Er 



ur 


iB/us 


( 2 . 2 . 2 ) 


where Ax is a step size for the scan. It should be noted that this attenuation includes not only 
1) so-called viscous loss due to interaction with thermally excited elastic waves, 2) scattering 
by crystalline defects, impurities, and surface roughness, and 3) air loading but also the losses 
due to diffraction and the leaky nature of the wave on GaAs [2.2.16]. The effects of loss due 
to interaction with electrons should be negligible because of the high resistivity of the 
substrate. 
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Figure 2.2.3 Propagation loss in the frequency range of 230-290 MHz. 
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The propagation loss shown in Fig. 2.2.3 was observed in the frequency range of 
230-290 MHz, by comparing several transverse scans on the metalized surface with the 
distance of about 200X and averaging the measured values. In fact, the diffraction was 
relatively large compared to that of bare GaAs due to the characteristic of the ZnO film to be 
discussed in the next section and the narrow aperture of the IDT (W=8.3 or lOX). With 
greater diffraction, more SAW energy can be lost to bulk shear since SAW propagation on 
this cut and propagation direction of GaAs is a leaky SAW mode. Employing an IDT with a 
wider aperture may reduce the value of the propagation loss in the measurement. 

Nevertheless, the narrow aperture was used in order to obtain enough diffraction for the 
velocity surface measurements to be discussed in the next section. 

In the laser probe system, the output signal of the laser probe, which is detected with a 
photo-diode, is filtered into dc and rf components. The dc component signifies the reflectivity 
of the sample under test, and the rf component is related to the normal component of the 
SAW particle displacement. While only the phase information of the rf component is required 
for the velocity measurement, its magnitude is used for the attenuation measurement. Since 
the magnitude of the rf component is proportional to the dc component, it is normalized by 
the dc value The accuracy of the magnitude measurement depends on, among other 
things, the flatness of the surface. However, the sample is not perfectly flat, resulting in a 
slight defocus of the laser probe spot and a variation in thus, precise measurement of 
attenuation is very difficult. Furthermore, due to the diffraction loss of the IDT, the maximum 
deviation of the attenuation measurement was up to a few tenths of a dB/ps depending on the 
scan distance. 
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The measured data shows that the 1.6 dc triode sputtered film with the passivation 
layer had the lowest loss and the 4 pm rf magnetron sputtered film without the passivation 
layer the highest attenuation, which is comparable to what has previously been reported 
[2 2.1],[2.2.17]. The passivation layer had a noticeable effect on both sputtered 
films-especially for the rf magnetron sputtered ones. The reduction m the attenuation is 
perhaps due to the fact that the passivation layer provides a better platform for the film 
growth by compensating for stress differences between the film and the substrate [2. 2. 2], In 
spite of the diffraction loss, it should be noted that the lowest value of attenuation for the 
films is not severely worse than that of bare GaAs substrate, which has the values of 0.323 

and 0.567 dB/ps, for 200 and 300 MHz, respectively [2.2.18]. 

The ZnO/GaAs structures may not have suitable temperature stability required for 
some applications due to positive temperature coefficient delays (TCDs) of both bulk ZnO 
(35 ppm/°C) and GaAs (50 ppm/°C). However, since the passivation layer of SiO^ 
(-70ppm/°C) has a strong negative TCD, improvement of the temperature stability may be 
achieved by adjusting the passivation layer thickness [2.2.2]. An experimental compensation 
of TCD of GaAs was reported by using a Au/SiO, film [2.2.19]. No measurements of the 
temperature characteristics have been attempted in this research. 


2.2.4 Velocity surface 

2.2.4. 1 Knife-edge laser probe measurements 

The substrates investigated for the velocity surface were 1.6, 2.8. and 4.0pm rf 
magnetron sputtered ZnO films with 0.1pm thick SiO, passivation layer and 1.6pm dc triode 
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Figure 2.2.4 Beam profiles are separated by 200X when X is 12 pm. The solid line is the near 
the IDT and the dotted line is for 200X ap^. (a) 1.6 ppi DC triode sputtered ZnO, (b) bare 
GaAs. 
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sputtered ones with and without 0.2^m thick SijN^ passivation when X was 12 pm. The 
typical transverse scan data measured on the free surface of 1.6 pm dc triode sputtered 
ZnO/GaAs and bare GaAs are shown in Figs. 2.2.4(a) and (b), respectively, where the SAW 
propagation is along the <1 10> direction of GaAs. The two transverse scans were separated 
by a distance of 200X. Note that the acoustic diffraction with the ZnO film shown in Fig. 
2.2.4(a) is much larger than that of bare GaAs, which clearly shows an focussing beam profile 
in Fig. 2.2.4(b). Thus, the velocity surfaces of these structures are investigated theoretically 

and experimentally to explain the difference in the diffraction. 

The theoretical velocity curves of the ZnO film with a 0.1 pm Si02 passivation layer 
are shown in Fig. 2.2.5(a) as a function of the ZnO film thickness normalized to X. The 
propagation directions are centered about the <110> axis. Note that the SAW velocity curves 
for the ZnO film thickness of 0.13, 0.23, and 0.33 X on GaAs are concave upward, and this 
is in contrast with that of the concave downward curvature for bare GaAs. The concave 
downward and upward curvature infer respectively focussing and non-focussing behavior, 

and support the diffraction shown in Fig. 2.2.4. 

The effect of the passivation layer with 0.1pm thick SiOj or 0.2pm thick S13N4 is 

shown along with non-passivated structures in Fig. 2.2.5(b) for the ZnO thickness of 1.6 pm. 
Since the passivation layer is treated as an isotropic media like the c-oriented ZnO film, its 
presence is not expected to affect the shape of the curve, but shifts the entire curve slightly 
toward faster values. Considering the fast SAW velocities of single SiO^ and Si3N4, one can 
easily expect that the curve would shift upward as the thickness of the passivation layer 
increases as long as the energy is trapped at the surface. 
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Figure 2.2.5 Velocity surfaces by Xaguerrc technique.* Values indicate thicknesses of ZnO 
normalized by h= 12 ^im. (a) ZnO film thickness with 0.1 iim SiOj layer, (b) 0.13X thick ZnO. 
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Figure 2.2.6 Experimental velocity surfaces, ZnO films'notmalizcd by 12 pm. (a) Effect 
of ZnO thickness and 0.1 pm SiOj layer, (b) Effect of passivation for 0.13X (1.6 pm) ZnO. 
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Applying angular spectrum of plane waves (ASPW) theory to the transverse scans 
produced velocity surfaces shown in Fig. 2.2.6. The usage of the laser probe system for the 
ASPW theory has been described in detail in Hunt and Hunsinger [2.2.20]. Note that there is 
close agreement between Figs. 2.2.5 and 2.2.6. 

The analysis of the velocity curve in terms of an anisotropy parameter gives an 
immediate insight into the diffraction of a propagating SAW. Using a parabolic 
approximation, the velocity of the pure mode direction (<110> of GaAs) can be expressed as 

v(0)=v„(l - bd^) (2.2.3) 

where v„ is the velocity at 0=0, 6 is an angle off the pure mode direction in radian, and b is 
the anisotropy parameter. The beam is spreading for b<0 and focussing for b>0 comparing to 
the isotropic case, b=0 [2.2.21]. As shown in Fig. 2.2.5, the parabolic approximation for the 
0.33A, thick ZnO film is not appropriate since the curvature is almost flat (i.e. isotropic) near 
the zero angle. Therefore, the values of b for the 0.13 and 0.23X thickness were obtained 
from a least-square's curve fitting and are listed in Table 2.2.3 for the comparison between 
the theory and the experiment. 

The inversion of the velocity curves was investigated numerically for the other 
materials. We arbitrarily chose GaAs, GaP, InAs, and InP since they are substrates with cubic 
symmetry, and AIN and ZnO as the piezoelectric film both of which have hexagonal 
symmetry. Material constants were available from Auld [2.2.22]. The inversion was observed 
for all of the combinations; e.g., it occurred at 0.05A. thick ZnO or 0.03 A. thick AIN both on 
GaAs. The curves became flat as the film thickness was increased due to the isotropic 
characteristic of the films. In addition, SAW modes exist only for relatively thin films of AIN 
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due to its much faster velocity (5608 m/s) than those of the substrates (1990-3540 m/s); ..g., 
for 0.1 a or less thickness of AIN on GaAs. For comparison, the velocity curve was also 
calculated for an isotropic film such as fused silica (3411 m/s) on GaAs, but in this case the 
inversion did not occur. The calculation results suggest that both the hexagonal films (AIN 
and ZnO) cause the inversion of the velocity surface in a cubic material for the cases 
mentioned above. This indicates that it is possible to control the diffraction of the substrate 
with proper selection of film thickness. 

2.2.4.2 Line-focus-beam scanning acoustic microscope measurements 

Again, Fig. 2.2.5(a) indicates that the velocity surface inversion occurs about the 
<110> direction when the c-oriented ZnO film is deposited on a {001)-cut GaAs wafer. Our 
theoretical calculations also reveal that the reverse situation appears if the propagation is along 
another pure mode axis <100> direction. In order to evaluate this predicUon and allow the 
measurement of acoustic properties along non-piezoelectric directions, velocity measurements 
using a Line-Focus-Beam scanning acoustic microscope (LFBSAM) 3 operating at 225 MHz 
were performed. The technique to deduce the velocity is the so-called V(z) analysis which has 
been reported in detail in [2.2.23]. Due to the requirement of the water coupling between the 
SAM lens and the sample, the velocity measured is associated with leaky waves for which the 
waves leak their energy into the water. Thus, the measured velocity due to the water loading 
is slightly different from the one measured by the laser probe. The mass loading of the water 
is included in the V(z) analysis. The LFBSAM can also measure both leaky SAW and 
pseudo SAW [2.2.23] which can exist in a [001}-cut cubic crystal. The advantage of this 
method is the ability to measure the velocity profile in 360 degrees. However, the spatial 
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Figure 2.2.7 Velocity surfaces measured by LFBSAM at 225 MHz. (a) Effect of ZnO film 
thickness with a 0.1 [im SiOj layer, (b) Effect of passivation for 1.6|im ZnO film. 
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resolution is limited (in this case an area of 300 |*m by 1000 ^m) because of the defocussing 
requirements [2.2.231. For the bare [001 )-cut GaAs, both leaky SAW and pseudo SAW 

profiles were reported [2.2.24]. 

Fig. 2.2.7(a) shows the measured velocity surfaces for the 0 (bare). 1.6, 2.8, and 4.0 
pm ZnO nim on (OOl)Ktut GaAs. It is seen that the ZnO ftlms change the focussing behavior 
about the <1 10> direction into the non-focussing one and the reveree inversions appear about 
the <100> direction. These results support our theoretical results and our laser probe 
measurements. Effects of the passivation layer (0.1 pm StO, or 0.2 pm SijN.) for 1.6 pm 
ZnO and |00I [-cut GaAs on the velocity surface are provided in Fig. 2.2.7(b), and they also 
agree with the results obtained using the laser probe. Note that Figs. 2.2.7(a)(b) arc solely 
used to demonstrate the leaky SAW behavior about the <1 10> and <100> directions, the 
details of their behavior farther away from these angles are not shown. 

Experimental velocity surfaces with and without the metallization pad for small 
deviations about the <110> direction are given in Fig. 2.2.8. Because the velocity 
measurement accuracy performed by LFBSAM is better than «.02%, in principle these data 
can be used to deduce the value of K’. Recently, the V(z) analysis of LFBSAM measurements 
for a multilayered anisotropic structure has been reported [2.2.25). However, due to a lack of 
understanding of the dielectric loading of the piezoelectric effect by the water, measurements 
of K’ can not be carried out at present. On the other hand, the water is distilled and 
de-ionized to have a resistivity of 10’«-cm. The conductivity effect of the water is mingled 
with its dielectric effect because the characteristic electromagnetic impedance is equal to the 
root square of the permeability over the permittivity. 
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Figure 2.2.8 Velocity surfaces measured by LFBSAM for small deviations about <11^ 
Effect of ZnO thickness with 0.1 lim SiO^ layer, (b) Effect of passivation for 1.6pm ZnO 


(a) 
film. 


65 




2.2.5 Reflection properties of metallic gratings 
An exact theoretical analysis of the reflection properties is extremely difficult. The 
reflectivity and velocity shift of metal gratings on a single crystal can be analyzed by using 
Datta and Hunsinger’s technique [2.2.3] if the metal thickness is small compared to the SAW 
wavelength. Considering that the metallic gratings effect the surface, Datta and Hunsinger s 
theory developed for a single material can be directly applied to a multilayered substrate by 
obtaining effective material constants of a single material equivalent to a multilayered 
substrate. The effective constants c,ff may be obtained by integrating multilayered material 
constants weighted with SAW energy distribution and normalized by the energy as follows, 

fayc(y)f(y) (2,2,4) 

f<fyP(y) 

where c(y) is a constant of the multilayered material tensor and P(y) is the energy 
distribution as a function of the depth y. Since the substrate is piezoelectric, the mechanical, 
electro-mechanical, and electrical energy distribution are utilized for stiffness, piezoelectric, 
and permittivity tensor, respectively. The detailed expression for the energy distribution is 
available from the Laguerre polynomial technique [2.2.7] and is not reproduced here. The 
validity of using the effective constants may be checked by comparing between the SAW 
velocity of ZnO/GaAs and that of the equivalent single material as is shown in Fig. 2.2.9. 

The maximum difference between both velocities is 0.009. 

The reflectivity of gratings can be divided into two components: the piezoelectric and 
the mechanical. The former depends on the piezoelectric coupling constant K^ and the latter 
depends on the thickness h of the metal in wavelength X. The reflectivity r of a metal strip 
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Figure 2.2.9 A comparison between SAW velocity of niultilayer (ZnO/GaAs) and that of 
effective single material. 
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with a metallization ratio n at a frequency f can be given by 




Sin 


i\nf 


(2.2.5) 


where f, is the center frequency, P. is the piezoelectric scattering coefficient, and F, is the 
first order mechanical scattering coefficient. P. is the same for all substrate-electrode 
combinations and depends only on the metallization ratio. On the other hand, F. depends only 
on the substrate and electrode materials and is independent of the metallization ratio. More 
detailed explanation about these parameters can be found in Hunsinger [2.2.3]. 

The value of P. at f„ is given as a funcuon of T1 by using the Legendre function 
(2.2.3). It is to be -0.75 when the value of q is 0.5. Fig. 2.2.10 shows the calculated values 
of K" and F. for an aluminum grating with 1/4 X width and 1/4 X spacing (the same 
geometry for solid finger interdigital transducers (IDTs)) on ZnO/0.25 pm SijNyOaAs as a 
function of ZnO film thickness. As shown in Fig. 2.2.10, the value of K= increases from 0.07 
% (GaAs) to 1.0 % (ZnO). and that of F. decreases from 0.9 (GaAs) to 0.2 (ZnO) as the film 
Ulickness increases. Thus, it is possible for the piezoelectric scattering and the mechanical 
scattering to cancel each other due to their opposite signs. This zero reflectivity can not only 
reduce triple transit interference (TTI) of IDTs for some filter applications but also allows 
the use of a single finger IDT instead of a split finger IDT, thus increasing the frequency 

range of easily fabricated devices. 

In order to verify this remarkable theoretical expectation, a set of characterization 
devices shown in Fig. 2.2.1 1 was designed to measure the effect of the gratings. The ZnO 
films were grown by dc sputtering method with thicknesses of 0.8 and 1.6 pm with a 0.25 
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Figure 2.2.10 Piezoelectric coupling and first order ihecbanical scattering coefficient of 
al uminum gratings on ZnO/0.25 fim Si 3 N 4 /GaAs. 
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First Order Mechanical Scattering 







(im CVD grown SijN 4 passivation layer. The probe IDTs are 100 "K wide where "k is 16 pm 
with a center frequency around 180 MHz. Each of the probe EDTs is 20 k long, and split 
finger electrodes are employed in order to reduce internal reflections. Both shorted and open 
gratings with the width of 108 ^ are investigated, and each of the gratings consists of 200 
strips with 1/4 X, width and 1/4 k spacing. The distance from the center of the front two 
IDT electrodes with same polarity to the edge of the shorted grating is designed to be " 
(l/2n+l/4) k and l/2nX, where n is an integer, to measure the sign of the reflectivity since 
the phase construction of a standing wave occurs at the distance of (l/2n+l/4) k for a 
negative sign of reflectivity and $l/2n \lambda$ for a positive one [2.2.26]. 

The properties of open gratings were also investigated along with the shorted gratings. 
Unlike the shorted gratings, the SAW regeneration from the open grating should be 
considered in the analysis. At the center frequency, the reflectivity, r„, of an open metal strip 
with the metallization ratio of 0.5 is given as [2.2.27] 

r^=r + jOMK^ (2.2.6) 

Due to this regeneration, the reflectivity of the open gratings always retains a positive 
value with a larger magnitude than that of the shorted gratings; thus, the open gratings may be 
suitable for reflectors or resonators, for which a strong reflectivity is desirable. 

The identical IDT probes with free and metalized surface with a 120 A, wide and 200 
k long metal pad instead of the gratings were also fabricated to calibrate the measured data. 
Moreover, from these calibration sets, the SAW properties of velocity, K^, and attenuation for 
the ZnO/GaAs substrates were also characterized. Three different thicknesses of the 
aluminum gratings were deposited with an e-beam evaporator. The measurements using a 
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surface profile measurement apparatus showed they were 400, 800, and 1200 A, which 
corresponded to 0.25, 0.50, and 0.75 % of X, respectively. The measurements were performed 

with a HP8753C network analyzer with a time-gating feature. 

The velocities for the 0.8 pm and 1.6 pm ZnO thicknesses were 2928 m/sec and 2912 

m/sec, respectively, from the center frequency measurements. The value of was also 
obtained by comparing the difference of the time delay from the IDT 1 and 2 to the IDT 1 
and 3 between the free surface and the metalized surface. It was 0.73 % and 1.1 %, 
respectively, but their accuracies may be degraded by the inaccurate value for the free surface 
velocity obtained from the center frequency measurements due to the velocity shift under the 
IDTs. More accurate measurement using the laser probe or the LFBSAM was not attempted 
for these substrates, but the measurements for the other substrates are already discussed in 
the previous section. The propagation loss was also measured by comparing insertion loss of 
the IDT 1 and 2 pair and that of the IDT 1 and 3 pair with a similar manner described in 
Melloch and Wagers [2.2.18]. Note that the loss due to diffraction is not considered in this 
technique. The loss of the 0.8 pm and 1.6 pm thickness was obtained to be 0.7 dB/ps and 0.8 
dB/ps on the metal pad, and 2.6 dB/ps and 3.6 dB/ps on the free surface, respectively. This 
difference between the free surface and the metal pad is due to the large diffraction of the 
ZnO/GaAs. The 120 X wide metal pad must act as a waveguide to confine the acoustic 
energy and reduce the insertion loss, resulting in the reduction of the propagation loss. Note 
that the results coincide well with the measurements by the laser probe, which are discussed 

in the previous section. 
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Figure 2.2.12 Typical frequency response of two shorted gratings with the same metal thickness 
showing phase construction and destruction. ' 


73 




Fig. 2.2.12 shows typical unmatched frequency response between the IDT 1 and 2 pair 
for the two shorted gratings (as mentioned before, the difference of their distance from the 
IDT 1 is 1/4 X.) with the same metal thickness. Note that the phase construction and 
destruction are clearly shown about f,. The reflectivity of the grating was obtained using the 
time gating technique described in Wright [2.2.27], Comparing the response of the IDT 1 and 
2 pair with the grating to that without gives a total reflectivity of the gratings assuming the 
distance between the IDT 1 and the grating is exactly correct for the total phase construction 
and destruction. Note that the propagation loss is negligible. The reflectivity of each electrode 
of the grating is given by [2.2.26] 

■ ■ tanh~^(/Q (2.2.7) 

' ' N 

where R is the total reflectivity and N is the number of the electrode. 

Figs. 2.2.13 (a) and (b) show the measured data and theoretical calculation for the 
reflectivity of an electrode in the shorted and the open grating at f„ as a function of the 
metal thickness normalized in X. Prior work concerning gratings on GaAs can be found in the 
literature[2.2.28]. The value of was assumed to be 70 % of the theoretical value based 
upon the measurements presented in the previous section. From the theoretical data, it is seen 
that for the shorted gratings the positive reflectivity on bare GaAs decreases, and eventually it 
becomes negative. The change in the reflectivity of the open gratings with film thickness is 
not so drastic as for the shorted gratings. On the other hand, the reflectivity on bare GaAs at 
f„ is given by 1.1 h/X, and its value has a marginal difference between the shorted and the 
open gratings due to the very weakly piezoelectric regeneration. From the measurements, it 
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Figure 2.2.13 Reflectivity of an aluminum electrode witb l/4X width and spacing 1/4X, at 
where the metal thickness in X. is represented by h: (a) shorted grating, (b) open grating. 
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should be noted ,ha, .he sign of .he reflec.ivi.y changes depending upon .he metal thickness 
and aero reflec.ivi.y occurs a. ,he metal duckness of 0.50 % k for 0.1 X (1.6 pm) dtickness 
of the ZnO nim for the shorted grating while the reflectivity for the open grating is always 
positive. The measured data shows relatively good agreement wifl. the theoretical prediction. 
The difference between the theory and the measurements may be due to the inaccurate 
effective material constants as well as the variation in the uniformity of the film and the 
metal thickness over the entire surface. 

The measurement of velocity shift of the grating was attempted from the shift of the 
center freqnency of the grating only response, ^ The shift of f., depending upon the metal 
thickness was not noticeable for both ZnO films, but the ftequency shift between the shorted 
and the open grating was observed with the value of 0.05T1.25 MHz. It is known that the 
regeneration of SAW from the open grating causes an increase m the velocity by 0.25 K at 
„=0.5 [2.2.27]; thus, the measured values are comparable to the calculated ones of 0.18 MHz 
and 0.25 MHz for 0.8 pm and 1.6 pm thickness of the ZnO film, respectively. 


2.2.6 Summary of Thin Film Work 

Theoretical and experimental data have been provided for the design of ZnO/GaAs 
based SAW devices, which may have the potential application for the monolithic integraUon 
of SAW and electronic devices. An analyttcal tool for SAW properties of a piezoelectric thin 
film was developed in this study, and its accuracy was 
of ZnO/GaAs substrates. 


verified with the experimental results 
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The SAW properties of both rf magnetron and dc triode sputtered ZnO films on 
{001) -cut <110>-propagating GaAs have been measured and compared to theoretical 
calculations. The passivation layer improves the ZnO film quality appreciably for both 
sputtering methods by reducing the propagation loss 0.5- 1.3 dB/ps at 240 MHz depending 
upon the ZnO film thickness. The best quality of film observed has been 1.6 pm dc triode 
sputtered ZnO film with 2000A SijN 4 , for which the propagation loss was 0.9- 1.0 dB/p"s at 
240-280 MHz. Its value of would be larger than that of the film without the passivation 
layer, 0.65-0.75%. Compared with bare GaAs, this particular film substrate has approximately 
0.6 % lower SAW velocity, at least 10 times larger K\ and 0.6 dB/ps more propagation loss. 

The velocity surfaces of the ZnO/GaAs substrates have been measured using a 
knife-edge laser probe and a LFBSAM and have been theoretically predicted. It is found that 
for certain thicknesses of ZnO film the velocity surfaces about the <1 10> direction are 
concave upward in contrast with the concave downward curvature of the bare GaAs. This 
causes the SAW to diffract as soon as it leaves the transducer. Theoretically, c-axis oriented 
hexagonal films such as AIN and ZnO can be used to control the diffraction of a cubic 
substrate such as GaAs, GaP, InAs, and InP by varying the film thickness. 

The reflectivity of shorted and open metal gratings has been analyzed and measured. 
The analysis showed that zero reflectivity could be obtained by a combination of the ZnO 
film and the metal thickness and the metallization ratio. Experimentally, zero reflectivity at 
the center fi-equency was observed at 0.5 % X thickness of a shorted aluminum grating with 
1/4 X width and 1/4 X. spacing on 1.6 pm thick ZnO film with 0.25 pm thick SijN 4 
passivation layer over GaAs. 
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TABLE 2.2.1: ZnO Film Deposition Parameters 


Sputtering method 

RF magnetron 

DC Triode 

Target 

ceramic 6.5 inch dia. 

ceramic 4 inch dia. 

Target to substrate 

1.5 inches 

2.5 inches 

Background pressure 

8 mTorr 

3 mTorr 

Substrate temperature 

350° C 

250° C 

Gas 

82% Ar 18% O 2 

90% Ar 10% O 2 

Power 

400W 

150W 

Rate 

4.6 pm/hour 

1.4 pm/hour 
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TABLE 2.2.2: Number of Finger Pairs N is 50.5 for all IDTs 


No. i 

Type 

Nominal 

Measured 

Aperture 

Aperture 



wavelength 

center 

(Wj) 

(W,) 



) [pm] 

freq.(fj 






[MHz] 

[pm] 

- - 

l\ ] 

1 

split 

16 

175.0-182.2 

80 

5 

2 

solid 

12 

230.2-240.6 

100 

8.3 

3 

solid 

10 

273.5-287.4 

100 

10 

4 

solid 

8 

339.0-356.8 

80 

10 
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TABLE 2.2.3: Anisotropic Factor of 0.13 and 0.232. Thickness of ZnO Film. The Values 
Are Obtained by Curve Fttting With a Uast-Snuare Method. The Values in the Parentheses 

Are for the Short Surface. 



0; None. 1: 0.1 ^lm SiO^ , 2: 0.2^m Si 3 N 4 
* Reference [2.2.11] 
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2.3 Acoustoelectric Amplifiers 


In the previous report, the concept of acoustoelectric amplification was introduced, an 
analytical model was briefly presented and the motivation for the research in this area was 
developed. In this section we present some theoretical predictions for bilevel structures using 
the model of Henaff et. al [2.3.1]. We also discuss the segmented structure which will be 
fabricated and tested over the first quarter of 1995. Finally the application of these devices to 
the HACT imager is discussed including size and power consumption considerations. 

2.3.1 Bilevel Model 

The simple model of [2.3.1 ] was adopted from bulk acoustic (BAW) wave amplifiers 
[2.3.2] and extended with great success to surface acoustics wave (SAW) devices. Due to its 
relative simplicity and accuracy this model found widespread use in modelling devices 
consisting of a semiconductor layer in close proximity or deposited on a piezoelectric 
substrate. It is also valid for a piezoelectric film deposited on a semiconductor or a bulk 
piezoelectric semiconductor. The details of the derivation may be found in [2.3.1]. 

Figure 2.3.1 shows schematically a amplifier structure on bulk n-type GaAs. The wave 
is launched and received by the interdigital transducers (IDTs) at the ends of the device and 
the surface wave is amplified due to the interaction with the drifting carriers accelerated by 
the voltage applied at the ohmic contacts. Figure 2.3.2 shows a comparison between the 
bilevel model and the experimental measurements on a 10'^ O 2 doped GaAs substrate at a 
synchronous frequency of 176 MHz [2.3.3]. Data points were taken from Figure 3 of [2.3.3] 
and smoothed using the SPLINE function in MATLAB to produce a curve for the comparison 
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Figure 2.3.1 Schematic diagram of the amplifier structure on bulk GaAs. 
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Figure 2.3.2 Comparison of the bilevel model and experimental measurements of GaAs 
substrate. 
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with the bilevel model. The y axis is the gain measured over a 1 cm separation of the ohmic 
contacts and the x axis is the voltage applied between the contacts. The model gives a good 
first order approximation to the experimental device and serves as a good tool to make 
predictions for our first experimental designs. The model deviates from the experiment at high 
voltage but more accurate models will be developed as required throughout the project. 

2.3.2 Segmented Amplifier 

The device of Figure 2.2.1 had several undesirable characteristics, the most prominent 
being the large voltage required across the ohmic contacts to amplify the acoustic signal.. This 
high voltage problem was resolved by using a segmented structure [2.3.4], The amplifier is 
segmented into many smaller sections over which a much smaller voltage is applied. The 
length of the segment is defined as LI and the segment separation as L2 as illustrated in 
Figure 2.3.3. The gain of the structure is the difference between the sum of the gain in the LI 
sections and the loss in the L2 sections. Such structures were demonstrated successfully on 
CdS with Ll= 500 urn and L2 = 250 urn [2.3.5], The voltage for maximum gain for these 
devices occurred around 120 V. 

Our approach is to design devices with the shortest possible segment length to give the 
minimum power dissipation. Devices have been designed with LI on the order one acoustic 
wavelength to enable the device to be operated with drift voltages less than 5 V. The bilevel 
model has been extended to predict the gain of the segmented structure. As a comparison to 
the high voltage device, the curve of Figure 2.3.4 shows the predicted performance of a 
segmented structure on 10'^ doped bulk GaAs with LI =100 urn and L2 = 20 urn at a 
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Figure 2.3.3 Proposed amplifier design using smaller segments. 
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synchronous frequency of 176 MHz. In order to achieve low voltage operation, some tradeoffs 
in size and power consumption must be made. As shown in Figure 2.3.4, the gain per unit 
length is smaller for the segmented structure than the high voltage device due to the loss 
incurred in the L2 sections. The power consumption may also be higher depending on the 
length of L2. These tradeoffs will be examined in more detail as the project progresses. 

Our devices will be fabricated on substrates consisting of a n-type epilayer on top of 
intrinsic GaAs. The device operation is the same as for bulk material and the bilevel model 
has facilitated the design of the first experimental pass. The layer stackup will be consistent 
with that used in the HACT imager so that the amplifier may be integrated into the imager as 
will be discussed in a later section. The structures intended for use in the HACT will consist 
of single LI gain sections placed throughout the imager. The segmented structure serves as an 
experimental device to develop an accurate model and is also intended for signal processing 
applications as a stand-alone device or integrated into a IC with other functions. 

2.3.3 ZnO Overlay 

The advantage of using a GaAs substrate for this device is the high carrier mobility. In 
the region where the carrier velocity is less than the saturation velocity, the carrier velocity is 
equal to the product of the mobility and the applied electric field. The gain of the device is 
proportional to the carrier velocity and so is proportional to the applied electric field (ie. 
applied voltage) times the mobility. Out of this it may be seen that the higher the mobility the 
lower the applied required to obtain the same gain. 
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The gain is also proportional to the electromechanical coupling coefficient, k^. The 
disadvantage to using GaAs is the relatively low value of k^ compared to substrates such as 
LiNbOj. However this problem may be overcome by using a ZnO overlay to increase the k^ 
by at least an order of magnitude. Figure 2.3.5 shows the gain vs voltage curve for a ZnO 
overlay with an estimated k' = 0.5% on the device modelled in Figure 2.3.4. We observe a 
peak gain improvement of approximately 35 dB with the addition of the ZnO film. 


2.3.4 Application to HACT Imager 

As mentioned in the previous report, the motivation behind this project is to enable the 
insertion of small gain blocks in between the pixels of the HACT imager to maintain a 
constant SAW amplitude across all of the pixels. Our plan is to use minimum length blocks to 
maintain minimum voltage operation. For example, if we use a single 8 um long gain block at 
360 MHz ( 1 wavelength) with the ZnO overlay the peak gain would be 0.07 dB and the 
power dissipation would be 0.036 mW. This by itself is not sufficient to overcome the SAW 
propagation the loss over the length of the imager. However, if for example we insert 100 of 
these miniature blocks along the propagation path of the SAW, then we would have 7 dB of 
gain at a power dissipation of 3.7 mW. The number of blocks will be optimized depending on 
the loss incurred over the SAW propagation path and the available real estate on the imager. 

2.3.5 Summary 

The bilevel model has been demonstrated for a bulk GaAs substrate. The model is 
equally applicable to multilayer piezoelectric structures. A segmented structure must be used 
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Figure 2.3.5 Gain versus voltage for ZnO overlay on GaAs. 
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in order .0 use this device with a practical voltage. Sorae gain and power consumption 
tradeoffs which result from the segmented structure will be examined later. The small gam 
sections will provide a method by which the SAW ampUtude may be maintained at a constant 
level across the length of the HACT imager resulting in an overall lower power consumption 

for the device. 
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2.4 Advances in ACT Device Modelling 

Two different approaches to the modelling of acoustic charge transport devices have 
been accomplished. The first is the optimization of the heterojunction acoustic charge 
transport device. Material layer thicknesses and doping concentrations were varied in an 
attempt to provide optimal characteristics. The second is modification of the existing code to 
simulate ZnO on silicon acoustic charge transport devices. 

2.4.1 GaAs Heterojunction ACT Optimization 
Previous reports detailed the computer model used for the simulation of acoustic 
charge transport devices and the design goals required for successful operation. One such 
device structure is displayed in Figure 2.4.1. These computer simulations of the layer 
structure displayed high charge carrying capacity. However, it was determined that due to the 
thin AlGaAs layer and the relatively deep location of the charge packet that non-destructive 
sensing contacts on the surface of the device would be highly inefficient. The image charge 
created by the signal electrons would appear in the substrate, as shown in Figure 2.4.2, 
instead of the metal layer. This is due to the distance from the electron packet surface being 
three times larger than the distance from the packet to the substrate electrons. In order to 
account for this short coming in the previous design a new design with a thick layer of 
AlGaAs was examined. Figure 2.4.3. 

Simulations of the new structure revealed the surface charge insensitivity of the 
previous design with essentially the same charge carrying capacity. 
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Figure 2.4.1 

state density. 


ilerial layer strocture of a HACT device insensitive to fluctuations in the surface 


92 



uouBJiueouoo uoJfO00 



Figure 2.4.2 Signal charge (electrons) in the channel and the associated image charge created 
in the GaAs substrate. 





FJonr.. 2 4 3 New material layer structure for aHACT device which is insensitive to the surface 
^r^Lnttation and will snppon image charge « the surface due to the thick AlGaAs layer. 
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2.4.2 ZnO on Silicon Acoustic Charge Transport 
Modifications of the device simulation code have enabled simulations of ZnO thin 
layers on silicon substrates. Since silicon is not piezoelectric, the acoustic wave is supported 
in the ZnO layer and the induced potential extends into the underlying silicon. However, 
since ZnO is highly resistive simulation of the layer within the device model would be very 
difficult due to the low carrier concentration. The modification of the code enables the 
simulation of the acoustic wave only in the ZnO, and both the acoustic potential and the 
carrier concentrations in the silicon substrate. Figure 2.4.4 displays the material structure, and 
Figure 2.3.5 shows a conduction band contour plot created in the silicon by the SAW 
potential. The SAW has created a potential well, very similar to that previously reported in 
GaAs, which is capable of transporting charge. 

2.4.3 Future Work on Device Simulation 

Future work in the area of acoustic charge simulation will require further optimization 
of the ZnO/Si system devices. In addition, future code modifications should make it possible 
to simulate moving charge packets by using a moving simulation domain solution. In this 
manner synchronis potential effects and injection phenomena maybe modelled leading to 
highly efficient device designs. 
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Figure 2.4.4 Material layout for ZnO overlay on a silicon substrate. 
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Figure 2.4.5 Conduction band equi-energy contours ia the silicon indicating a potential well 
which is capable of transporting charge. 
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2.5 Materials Growth 


In previous reporting periods materials growth development has focused on the 
development of an indium-free mounting technique for 2 inch substrates along with 
investigation of the spatial uniformity of dopant and host fluxes. In the present reporting 
period efforts have concentrated in three areas; substrate heater modification, system cleaning 
and temperature controller upgrading. 

2.5.1 Substrate Heater Modification 

Investigations into the substrate temperature uniformity showed that large temperature 
variations across the substrate surface were present. These temperature variations were due to 
the four filament design of the substrate heater in the Varian Gen II molecular beam epitaxy 
system. Even upon rotation and the incorporation of a pyrolitic boron nitride diffuser between 
the heater elements and the substrate rear surface, the effects of temperature vanations could 
not be eliminated. Therefore, a completely new substrate heater was designed, tested and 
incorporated into the existing substrate heater mount. This heater is based upon a broad area 
pyrolitic boron nitride encapsulated graphite heater element. Not only has this heater lead to a 
much more uniform substrate temperature, but is also mechanically more robust than the 
original Varian design. This design should eliminate substrate heater failure, which in the past 
has been the major reason for system down time. 
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2.5.2 System Super Cleaning 


The second major effort in the past six months has been the complete cleaning of the 
system. It was determined that the surface particulate density was higher than acceptable for the 
large area devices required for this program. The majority of these particulates were found to 
be generated in the growth chamber. To reduce the particulate density the entire system was 
opened and scraped out. This has resulted in an order of magnitude improvement in the system 
base pressure, along with a shortened bakeout time. A reduction in surface particulate density by 
a factor of four has be observed after the system cleaning. 

2.5.3 Gallium and Aluminum Controller Replacement 

Recent publications along with our previous investigations lead us to closely monitor our 
A1 and Ga effusion cell temperatures and fluxes. Shown in Figure 2.5. 1 is a typical plot of the 
gallium effusion cell temperature as a function of time with the original Eurotherm model 919 
temperature controller and model 932 SCR. One startling feature is that the temperature 
undergoes periodic variations with a period of about 400 seconds. This is similar to that which 
has been observed in other research laboratories. Not only are these p>eriodic variations in 
temperature observed, but long term variations in temperature of up to +/- 0.5 °C are also 
observed. 

To improve the temperature and therefore host element flux stability the older Eurotherm 
temperature controllers and SCRs were replaced on the aluminum and gallium cells. A plot of 
the typical time dependence of the of the gallium effusion cell using the new Eurotherm model 
818s temperature controller and model 832 SCR is shown in Figure 2.5.2. As can be seen the 
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periodic temperature disturbances have been essentially eliminated and the total variation in 
temperature has been reduced to +/- 0.1 °C. This has resulted in a very stable effusion cell flux 
as shown in Figure 2.5.3 where the time dependence of the aluminum cell flux, as monitored by 
an ionization gauge which has been rotated into the growth position, is shown over a period of 
1800 seconds. The total variation of the flux has be reduced to less that 2 %. 

Another effect leading to nonuniform doping and/or growth rate that has been investigated 
is the flux transients that are caused by changes in effusion thermal loading when the effusion 
cell shutters are opened to initiate growth. The aluminum flux as a function of time is shown in 
Figure 2.5.4, as the shutter is opened. As can be seen during the first 200 seconds of growth (~ 
450 A) the aluminum flux decreased by 14%. Not only are the flux transients observable upon 
the initiation of single layer growth but are also seen during superlattice growth. The gallium flux 
as a function of time during the growth of a GaAs/AlAs superlattice structure is show in Figure 
2.5.5. In this growth sequence the Ga shutter was open for 100 sec. and then closed for 100 sec. 
Upon each shutter opening in the structure a flux transient is observed. These flux transients lead 
not only to slight inaccuracies in the layers thickness, but also to small variation in the doping 
densities near the interfaces. At present feed-forward temperature control algorithms are being 
developed to eliminate these flux transients. 
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Figure 2.5.3 Stable effusion cell flux obtained using new controllers. 
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I Figure 2.5.4 Aluminum flux as a function of time after initially opening shutters, flux decreases 

I by 14%. 
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Figure 2.5.5 Gallium flux transients observed when the shutters are periodically opened and shut 
for 100 sec intervals. 
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3. Charge Transfer and Overflow Device 

3.1 Introduction 

Of the three devices that comprise the acoustic charge transport (ACT) imager, the charge 
transfer device (CTD) plays a major role in its operation, where the CTD collects and stores the 
electrons generated in the avalanche photodetector (APD), then transfers this accumulated charge 
to the ACT device for readout. The charge storage and overflow properties of the CTD directly 
affect image quality, while the transfer capabilities in part detemune the speed and fidelity of 

image readout. 

In this section the progress made on the CTD will be presented. An extensive review of 
charge storage device technology is presented, along with its application to the CTD. The 
electrical characteristics of the earlier vertical CTD (V-CTD) structures have been explained in 
light of other GaAs charge storage devices, as will be described to follow. This work on the V- 
CTD is also relevant to the development of the lateral CTD (L-CTD), since vertical charge 
confinement is necessary in either architecture. This section will conclude with future goals for 

CTD research and imager development. 

3.2 Review of Charge Storage Devices 

3.2.1 General Concepts for Charge Storage Device 
From a review of Si CCDs [3.2.1], two modes of charge storage are used in solid state 
imaging devices, the surface channel or enhancement approach and the buried channel or 
depletion approach. The charge storage operation in either mode is the same; they only differ in 
how the charge storage region is created. The features common to both modes are introduced 
first, followed by a description of the different charge storage regions in the surface channel and 
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buried channel devices. In the following discussions, the stored charge could be either electrons 
or holes, however, only electron storage will be described since the APD will operate with 
electron multiplication. 

In common to both storage modes, there exists a steady state concentration of electrons 
in the charge storage region or well of the device, where this electron concentration represents 
the full well capacity of the storage device. In this steady state, full well condition the device is 
not necessarily in equilibrium, i.e. it may be biased. Through the application of appropriate biases, 
a non-steady state, transient condition is created where the electrons in the storage region are 
removed to yield an empty well ready for charge storage. Electrons can be introduced into the 
storage well either by lateral injection or through optical generation within the well, until the 
electron concentration reaches its steady state value. Since the storage device is operating under 
transient conditions, the storage well will eventually return to the steady state, full well condition 
even if no electrons are purposely introduced. This degradation of the charge storage well occurs 
either by leakage of charge into the charge storage region, such as the reverse bias leakage across 
a p-n junction, or via spurious generation of electron - hole pairs, possibly aided by defects, 
within the well. These so-called dark current mechanisms determine the maximum storage time 
for charge in the device. 

Charge storage devices differ in how the potential well is created. In the surface channel 
or enhancement mode device, for example a Si MOS capacitor, charge is stored in an inversion 
layer at the interface between a semiconductor and insulator. This inversion layer, in turn, is 
surrounded by a depletion layer which acts to confine the electrons both vertically and laterally. 
The charge capacity of the surface channel device is determined by the bias applied to the gate 
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and by the capacitance of the dielectric. In a buried channel or depletion mode device, charge is 
stored in a potential well verUcally defined by two space charge regions. These space charge 
regions can be produced by p-n juncUons, an MOS structure, or a schottky barrier. An example 
of an depletion mode device is a GaAs MESFET structure, where the schottky gate contact has 
been biased to pinch-off the channel. Table 3.1 shows some general advanuges and disadvantages 

for these two charge storage modes. 


TmMP ^1: c:omparison ot Charge Storage Modes tor SoEd 
advantages 

Enhancement mode * largest possible charge storage 

capacity 

• simpler variation between Q and C 
compared to depletion mode 

• devices can be tested with no lateral 
isolation mechanism, such as a 
junction or mesa, because of 
depletion layer surrounding the 
inversion layer 

• doesn't require a good insulator - 
semiconductor interface, can use 
Schottky barrier or p - n junctions 

• fast charge transfer since low 
surface states 

• can use positive, negative, or both 
bias schemes 

• more radiation tolerant 


Depletion mode 


State Imaging Devices7 

DISADVANTAGE^ 
need the proper insulator - 
semiconductor interface, 
none easily implemented 
currently for GaAs 

slower charge transfer due 
to interface states 

other standard problems 
associated with MOS 
devices in Si due to 
presence of oxide 
leakage current often limits 
storage time rather than 
generation; especially true 
for GaAs schottky barriers 

smaller charge capacity 
compared to enhancement 
mode 

may require an extra 
contact to pinch-off or 
totally deplete structure 

need device isolation for 
test structures, probably 
mesa or schottky barriers 
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3.2.2 Charge Storage Devices in GaAs 

As noted in the previous section, there is no suitable insulator for GaAs that would allow 
enhancement mode charge storage devices based on a GaAs - dielectric interface, although 
progress is being made using novel structures grown by MBE [3.2.2]. Therefore, GaAs charge 
storage devices, such as CCDs, are based on the buried channel approach [3.2.3]. The standard 
GaAs CCD is based on a MESFET structure on a semi-insulating substrate, and is typically used 
for high speed signal processing and benefits from the easy integration of MESFET based 
circuitry with the CCD. 

An important consideration for the CTD is the storage time of the device. As discussed 
in the previous section, the mechanisms responsible for dark current degradation of a charge 
storage time are generation in the depleted well, which for GaAs should be a slow process 
because the generation rate is proportional to n„ and leakage currents. In such devices, leakage 
is the main well degradation mechanism, because of the small barrier height for schottky barriers 
on GaAs coupled with the small generation rate for GaAs. To assess the maximum leakage 
current tolerable by the CTD for HDTV applications, a simple calculation has been made 
following Fossum et al. [3.2.4], assuming that the storage time is leakage-limited, the maximum 
storage time necessary for the CTD is 5 ms (170 frames per second), and the well capacity is 10® 
e cm ^ a value consistent with current Si CCD-based HDTV imagers. This calculation yielded 
a maximum leakage current of 1.6 nAcm'^ This value should be compared with the lowest 
reported leakage current reported by Fossum et al. for 2DEG-CCDs of 30 g Acm'^ and a Rockwell 
result for a MESFET structure using an undoped AlGaAs layer of 0.01 - 0.1 nAcm l Although 
the Rockwell value appears encouraging, their results were for large area devices, where surface 
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effects are less significant. For the necessarily small area devices needed in the HDTV imager, 
such low leakage currents may not be possible. Therefore, an important issue that must be 
addressed for the CTD is development of low leakage current structures. 

GaAs CCDs are based on both standard MESFET and high speed HEMT (2DEG-CCD) 
structures [3.2.4]. Design constraints represent a trade-off between storage capacity and transfer 
efficiency, where increasing the storage capacity yields a large potential trough between the wells 
which captures signal charge during transfers. This problem is not as severe for Si CCDs because 
of the use of oxide passivation between the wells and overlapping gate structures. The minimum 
energy of this potential trough between wells increases with increasing channel doping and 
decreasing depth of the well from the surface, and decreases with decreasing well separation. 
Typical storage capacities for MESFET CCDs are around 1 - 5 x lO" cm , while 2DEG-CCDs 
have larger capacities around 1x10'^ cm ^ However, because the channel in a 2DEG-CCD is 
located very near the surface, the charge transfer efficiency suffers unless the wells are quite 
close (~ 1 Mni). This factor may not be important for application of 2DEGs to CTDs, since only 
one transfer is necessary between the CTD well and the HACT channel, and resistive gate 
technology could be used to remove the inter-well potential troughs [3.2.4]. A greater limitation 
to the use of 2DEGs for CTDs is the large leakage current of the schottky barrier gate due to the 
high electric field caused by the n^ AlGaAs layer [3.2.4], an issue that can be addressed through 
the use of p-i-n doped 2DEGs or improved schottky barrier technology. 

Charge storage devices have also been developed which are based on depletion mode 
storage using p-n junctions, where very long storage times have been observed [3.2.5, 3.2.6]. 
However, in these devices no free charge is transferred between storage wells, but rather only a 
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capacitance change due to the depleted charge is sensed, for use as digital memory devices. 
Therefore, these structures are not directly applicable to the CTD. 

Although a suitable insulator is not readily available for GaAs, enhancement mode devices 
have been attempted utilizing the heterojunction discontinuities between GaAs and AlGaAs, 
where the AlGaAs barrier is undoped. This structure is used in the heterostructure insulator gate 
FET (HIGFET) device. Unfortunately, because of the small energy barriers at the GaAs - AlGaAs 
interface (AE(.= 0.263 eV, AE^- = 0.132 eV for x = 0.3 AlGaAs) [3.2.7], the gate leakage currents 
are very large, in low milliampere range at room temperature [3.2.8], which would result in very 
small storage times. Kleine et al. [3.2.9] at Purdue have examined enhancement mode operation 
of an GaAs - AlGaAs structure, and only observed measurable electron inversion charge below 
70 K, where charge loss was limited by thermionic emission across the heterojunction conduction 
band discontinuity. 

One other GaAs storage approach that has been investigated is a floating gate device, 
where charge is stored in an GaAs - AlGaAs quantum well [3.2.10]. However, some question 
exist about the mechanism for charge storage in these devices, and that trapping effects may be 
responsible [3.2.1 1]. 

3.3 Results 

3.3.1 Comparison Between Simulations and Experimental Results for CTD-5 

The capacitance-voltage (CV) and transient capacitance (Ct) data from the V-CTD 
structure CTD-5 has been described based on simulations using Atlas II and HETMOD, a device 
simulator made available to Georgia Tech by IBM. Specifically, the storage capacity for CTD-5 
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has been calculated and aspects of the CV and Ct data explained. Simulation are currently being 
performed on CTD-8, with results to be presented in dte next report. This work has been 
hampered by the high failure rate of CTD-8 devices, as will be described in section 3.3.2. 

n,e equilibrium band diagrams for the two different V-CTD device structures are shown 
in Figures 3.1 (a) and (b), where we assumed originally that the charge storage well was formed 
by the AIGaAs layer and either a volume-doped (CTO-5) or delta-doped (CTD-8) p-type layer. 

To determine if charge storage could occur in the well, the effect of bias on the electron 
concentration in the well was simulated for CTD-5. For charge storage to occur, the free electrons 
must be removed from the well. The integrated electron density in the well was 2.7 x 10” cm” 
in equilibrium. However, no change in electron density, either an increase or decrease, was 
observable for readout (surface) biases from -1.0 to 1 JO V. Therefore, the well formed by the 
AIGaAs layer and the p barrier layer does not function as a charge storage well. The lack of 
charge modulation in this region results from all the applied bias being dropped across the two 
depletion layers, so that there is no potential variation and thus no charge modulation within this 

well. 

However, charge storage was observed in CTD-5 from the Ct measurements described in 
previous reports. The question is where in the structure is charge being stored? We can 
understand the charge storage mechanism operating in CTD-5 from an analogy to the GaAs 
MESFET depledon mode CCD, where an n-type channel is depleted by a suitable bias to pinch- 
off. In CTD-5, this channel is the p-type barrier region, with its barriers formed by the depletion 
layers of the two p-n junctions. During a Ct measurement, a bias is applied to the readout which 
depletes holes from the p-type region; when this bias is removed the p-type region fills with holes 
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Figure 3.3.1 (a) Equilibrium band diagram of a delta doped CTD structure (CTD-8). 
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Figure 3.3.1 (b) Equilibrium band diagram of a volume doped CTD structure (CTD-5). 
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via generation and leakage across the p-n junctions. Therefore, this structure acts as a hole 
storage device. Using HETMOD a pinch-off voltage and hole storage capacity for CTD-5 were 
calculated of 1.49 V and 1.15 x 10'^ cm ^ respectively. 

The storage lifetime (1/e value) observed in CTD-5 from Ct measurements at bias pulses 
slightly less than required for breakdown (± 4 V) were Tj = 3 ms. The lifetime is too short for 
the CTD, and is also small compared to literature values on the order of minutes, which indicates 
that either the material has a large number of defects which increases the generation rate in the 
p-type region or that leakage currents across the p-n junction barriers are large. From previous 
work we had identified that the leakage currents for CTD-5 were sizable and resulted from 
shunting of the front n-p junction due to excessive alloying from the AuGe - based ohmic 
contacts. Therefore, with the improved PdGe contact technology the charge storage time should 
increase. 

3.3.2 Device Characterization of CTD 

Investigation of the failure modes for V-CTDs using CV and conductance- voltage (GV) 
measurements have been performed as a function of device area, which isolates the effect of bulk 
phenomenon from mesa sidewall surface effects on device breakdown. For the volume-doped 
CTD-5, device failure resulted from ohmic contact diffusion into the junction region, as described 
previously. Therefore, an alternative ohmic contact technology for the V-CTD was developed 
based on PdGe metallization, which was used in all subsequent devices. For the delta-doped 
CTD-8, the major contributor to breakdown was a bulk-related effect, which we believe resulted 
from electric field induced failure of the Be delta-doped barrier layer due to movement of the 
dopant ion [3.3.1]. Future delta-doped structures will utilize low temperature. As rich growth 
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conditions to minimize Be diffusion, thereby improving the lifetime and operating characteristics 
of these devices. 

Initial investigation have begun on the effect of traps on CTD operation. The variation 
in capacitance with test frequency is being studied, where traps are not expected to contribute to 
the measured capacitance with increasing frequency, thus causing a decrease in capacitance. 

Storage lifetime measurements have recently been performed in liquid nitrogen (LN 2 ) at 
77 K, to study the effect of temperature on leakage currents and generation rates. Significant 
enhancements in the storage lifetime for CTD-5 have been observed, with a value of 68 s for a - 
1 V bias pulse at 77 K, compared to 2.9 ms for a - 4 V bias pulse at 300K. Future work will 
examine the dependence of storage lifetime on the bias pulse, pulse duration, and temperature, 
as a precursor to deep level transient spectroscopy setup and testing. 

3.3.3 Fast Storage Lifetime Testing System for CTD 

An approach for determining the very short charge storage lifetimes, of the order of 
microseconds, is to implement a fast transient testing system based on voltage sensing, rather than 
capacitance sensing, since the minimum capacitance measurement time is at best 10ms, compared 
to voltage measurements on a digital oscilloscope with times less than 1 ps. In such an approach 
a voltage proportional to the RC time constant for the CTD is measured after a pulse bias, and 
the dependence of device capacitance on time is deconvoluted from the time-dependent RC 
product. This approach actually determines a dc capacitance, as opposed to a standard capacitance 
meter which measures an ac capacitance. 
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Such a fast voltage transient system has been modeled at GTRI using the circuit 
simulation program PSPICE, to define a test setup where the CTD capacitance and a known 
resistance determine the RC time constant of the measurement setup, rather than any parasitic 
capacitances and resistances associated with the setup. This system has recently been setup, and 
is now being tested at GTRI. 


3.3.4 Summary and Future Work 

Significant progress has been made in our understanding of the operation of previous 
CTDs. This knowledge will allow us to more effectively design device structures using the 
HETMOD and Atlas II simulators for subsequent growth , fabrication, and testing. From the 
review of GaAs charge storage devices and the experimental results on CTD-5, an important need 
for long storage time CTDs is a low leakage current device structure, which can be realized by 
a high quality GaAs - insulator interface. Such a structure would benefit both depletion mode 
devices by reducing leakage currents, and also allow the fabrication of enhancement mode 
devices. Such a technology would also find application in many other areas of GaAs device 
technology. Work will also begin on construction of a deep level transient spectroscopy system 
for enhanced characterization of CTDs, APDs, and general MBE grown material. 
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4.0 Avalanche Photodiode Development 


4.1 Theory of Impact Ionization and Modeling of Avalanche Photodiodes 

One of the keys to developing an accurate model of an avalanche photodiode for device 
simulation is to formulate the impact ionization transition rate correctly. Most early theoretical 
studies of the impact ionization rate have limited validity since they neglect the complexity of 
the phonon scattering processes and the energy band structure of the host material. Numerical 
transport models based on the ensemble Monte Carlo method have been more recently employed 
[4.1. 1-4.1. 3] which include the full details of the band structure based on either a pseudopotential 
or k»p calculation. The impact ionization transition rate is typically formulated within these 
numerical transport models using a parameterized expression called the Keldysh formula [4.1 .4]. 
The Keldysh formula contains two parameters, p and E^,, which are very difficult to determine 
from first principles. Both of these parameters are selected such that the calculated ionization 
rate agrees well with experimental measurements. Often, more than one set of parameters can 
reproduce the same ionization rate in a given material. Subsequently, the nature of the ionization 
process cannot be uniquely determined from use of the Keldysh formula alone. In addition, the 
Keldysh formula is derived assuming several simplifying assumptions. These are that the 
transition is calculated assuming a direct semiconductor with parabolic bands for all of the 
carriers, no k-dependence is considered, and the overlap integrals are constant for all initial 
energies. 

There have been some recent attempts by several investigators to address the limitations 
imposed by the use of the Keldysh formula [4. 1.5-4. 1.9]. During the course of this project, we 
have developed a completely numerical study of the interband impact ionization transition rate 
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for both bulk silicon and GaAs [4.1.8, 4.1.9]. The ionization transition rate is evaluated using 
a wavevector dependent, semiclassical formulation of the transition rate. The transition rate is 
determined using Fermi’s golden rule from a two-body screened Coulomb interaction assuming 
energy and momentum conservation. The overlap integrals in the transition rate expression are 
determined through use of numerically generated wavefunctions determined via a k®p calculation. 
The full details of the calculation and the model are reported in references [4.1.8] and [4.1.9] 
which are enclosed within the appendix. It is found from these calculations that the transition rate 
is greater for electrons originating within the second conduction band than the first conduction 
band and that the threshold can be classified as relatively soft in both GaAs and silicon. It is 
further found that the threshold is somewhat harder in silicon than in GaAs. 

Though the numerical treatment of the interband impact ionization transition rate is useful 
it is highly computationally intensive. As such calculations of the transition rate typically require 
many cpu hours. In an attempt to reduce the computational demands necessary for a completely 
numerical formulation of the transition rate, we have investigated the usefulness of an alternative 
analytical formulation of the transition rate [4.1.10]. The new wavevector dep>endent analytical 
formulation was first developed by Quade et al. [4.1.1 1]. The details of their model have been 
sketched out in reference [4.1.10] which is enclosed within the appendix. It is found that the 
wavevector dependent model of Quade et al. produces very similar results to those obtained with 
the completely numerical model for some quantities. Specifically, both models predict that the 
effective threshold for impact ionization in both GaAs and silicon is quite soft, that the majority 
of ionization events originate from the second conduction band in both materials, and that the 
transition rate is wavevector dependent. Nevertheless, there are significant limitations to the 
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usage of the Quade formula for modeling interband impact ionization. These are that the Quade 
formula does not accurately reproduce the quantum yield data for bulk silicon, it requires 
determination of a new parameter, related physically to the overlap integrals of the Bloch states, 
and fails to account for any wavevector dependence of the overlap integrals themselves. As such 
the transition rate can be seriously overestimated at those poinu within the Brillouin zone for 
which nearly vertical, small change in wavevector, transistions occur. Based on these results, we 
have concluded that the Quade formula may not be truly useful in avalanche photodiode 
simulation and have retained the completely numerical technique for use in the neat generation 

of APD device simulators. 

In addition to our work on developing an improved theory of interband impact ionization, 
we have further worked on refining our APD device simulator as well as work on studying the 
performance of APD devices. Towards this end, we have begun work on tailoring our existing 
hydrodynamic simulator to study APD structures. This work has progressed along two related 
directions. First, we have fashioned the simulator for studying cylindrical geometries. This has 
been used to study a separate absorption, multiplication APD used by Bell Northern Research for 
lightwave communications systems. In addition, we are presently Investigating the inclusion of 
nonparabolicity effects into the hydrodynamic model to greatly enhance its accuracy for high 
energy devices. Two different formulations are being investigated, a simple nonparabolicity model 
and a power law model. Though the power law model offers somewhat better accuracy at high 
energy, its use is somewhat frustrated by the need to evaluate some of the moments numerically. 
We expect that we will soon incorporate the nonparabolicy technique into out existing 
hydrodynamic simulator enabling its usage for APD device simulation and development. 
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4.2 Characterization of Avalanche Photodiodes 


During this report period, a paper entitled "Gain Properties of Doped GaAs/AlGaAs 
Multiple Quantum Well Avalanche Photodiode Structures" was submitted and accepted for 
publication in the Journal of Applied Physics Letters. The paper presented experimental 
confirmation of a structure-induced carrier multiplication in the low-voltage regime due to 
interband impact ionization in the MQW APDs. This was contrasted with data from a 
conventional p-i-n APD where no such gain existed. The details of the experiments and the 
results were described in the previous report. 

Recent work has concentrated on performing several device simulations to help gain some 
insights into the physical phenomena that take place in p-i-n's, undoped, and doped MQW 
structures. Various simulations using the Gummel and Newton's methods were used to calculate 
electric field distribution, valence and conduction band potentials and carrier concentrations, as 
well as electron and hole mobilities throughout the structure. 

Figures 4.2.1(a) and 4.2.1(b) show the electric field profile and band potentials for a p-i-n 
structure. Both the p- and n- layers are heavily doped (1x10'* cm'*), while the i-layer is made up 
of intrinsic GaAs. Note the two peaks in the field plot which are the result of band bending at 
the edges due to bandgap narrowing caused by heavy doping. Figures 4.2.2(a) and 4.2.2(b) 
present the carrier concentrations and mobilities throughout the p-i-n structure. Electron and hole 
mobilities are lower in the doped GaAs regions due to ionized impurity scattering. 

Figures 4.2.3 and 4.2.4 show the results of the simulation when applied to an undoped 
GaAs/AlGaAs MQW structure with four built-in wells, while Figure 4.2.5 show the effect of 
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Figure 4.2.4 Simulated carrier concentrations and mobilities as a function of depth obtained for 
the MOW APD structure. 
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i Figure 4.2.5 Simulated electric field and band diagrams for a MQW structure with a 50 A wide 
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introducing two 50 A p-n doping layers in the center of every GaAs well. Note, how in Figures 
4.2.5(a) and 4.2.5(b), the doping in the wells creates peaks in the electric field and band bending 
in the potential wells which help increase the kinetic energy of the carriers and thus enhances 
impact ionization. This clearly illustrates the advantages of using doped MQW structures to 
greatly enhance the gain characteristics of such devices over that of undoped MQW APDs. The 
simulation results shown in Figures 4.2.6 and 4.2.7 show a quantitative measure of the 
relationship between the electron energy gain as a function of the doping density and the doping 
layer width. In Figure 4.2.6, the p and n doping layer widths were held constant at 50 A while 
the doping concentration was doubled (1.5xl0‘* to 3.0x10'*). The resulting gain in energy has 
increased by two fold from about 0.06 eV in Figure 4.2.6(a) to 0.12 eV in Figure 4.2.6(b). In 
Figure 4.2.7, the same simulation was conducted holding the doping density constant while 
increasing the p and n layer widths from 100 A to 150A. Once again, the energy band offset has 
doubled from about 0.45 eV in Figure 4.2.7(a) to 0.9 eV in Figure 4.2.7(b) resulting in a 
substantial increase in the kinetic energy of the electrons. 

The simulation shown in Figure 4.2.8 was conducted as an attempt to study the impact 
of doping imbalances on the field profile and the band structure. Such a problem occurs quite 
frequently during material growth and is not easily controlled. In this case, we chose the n 
doping in every GaAs well to be 20% higher than the adjacent doping in the p layer. The 
resulting electric field profile, as shown in Figure 4.2.9, gradually decreases. This results in 
partially depleted stages as shown in the band diagram in Figure 4.2.10. Partially depleted and 
fully undepleted stages in a MQW APD play a major role in limiting the gain characteristics of 
the device at low bias. Therefore, in order to fully maximize gain at low bias, it is very important 
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Figure 4.2.7 Energy gain for the p^n junction MQW structure when the thickness of the p-n 
junctions are increased to 100 and 150 A respectively. 
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Figure 4.2.8 Structure simulated to determine the effect of a doping imbalance in the p-n 
junction of the MQW structore. 
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Figure 4.2.9 Electric field profile obtained for the structure with the doping imbalance. 
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Figure 4.2.10 Band diagram of the 
stages of the device. 




to be able to deplete as many stages as possible which requires the need for precise doping 
control during MBE growth. As was previously reported, we were recently able to successfully 
fabricate several 10-period doped-well APDs that demonstrated full depletion below avalanche 
breakdown. However, when the number of built-in GaAs wells was increased to 25, the APDs 
displayed only partial depletion. This prompted us to try and determine the threshold limit beyond 
which further depletion ceases to occur. Therefore, we are currently in the process of growing 
and fabricating additional doped-well APDs where the number of stages in the avalanche region 
ranges from 10 to 25. This will help us determine the optimal width and periodicity of the 
depletion region which will maximize the gain at low bias in these structures. 

For future work, we will extend our simulation further to gain more insights into the 
properties of these structures at low and high bias when impact ionization becomes important. 
This will eventually lead to constructing better theoretical models of the gain, noise and 
bandwidth characteristics of the APD which can be used to analyze our existing experimental 
data and enable us to use such models as a basis for fabricating high performance APD 
structures. 
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* Ensemble Monte Carlo calculations of the electron interband impact ionization rate in bulk GaAs 
are presented using a wave-vector (k)-dependent formulation of the ionization transition^ rate. The 
transition rate is evaluated through use of numerically generated wavefunctions determined via a 
k*p calculation within the first two conduction bands at numerous points within a finely spaced 
three-dimensional grid in k space. The transition rate is determined to be greatest for states within 
the second conduction band. It is found that the interband impact ionization transition rate in bulk 
GaAs is best characterized as having an exceedingly “soft” threshold energy. As a consequence, the 
dead space, defined as the distance over which the ionization probability for a given carrier is 
assumed to be zero, is estimated to be much larger than that estimated using a “harder” threshold. 

These results have importance in the design of multiquantum-well avalanche photodiodes. 


I. INTRODUCTION 

The probability that an interband impact ionization event 
will occur depends upon the likelihood that a carrier will first 
drift to an energy equal to or greater than the band-gap en- 
ergy, often referred to as the threshold energy, and then suffer 
an impact ipnization event before phonon scattering pro- 
cesses relax it back to an energy below threshold. Most of 
the early analytical theories of interband impact ionization**^ 
concerned themselves with determining the probability that a 
carrier survives to threshold assuming that it would ionize 
immediately afterward. Such theories are said to assume a 
“hard threshold” for interband impact ionization. 

The hard threshold model was not limited to the early 
analytical theories. More advanced theories^'^ were formu- 
lated using a hard threshold model for interband impact ion- 
ization; however, more recent numerical studies^’^ as well as 
analytical theories^’* have called this assumption into ques- 
tion. Interestingly, in most of the numerical theories em- 
ployed in the past the calculated impact ionization rates often 
agree well with experimental results for some semiconductor 
materials using either a hard or “soft” threshold model.^ This 
is due to the fact that most numerical formulations rely on a 
parameterized expression, called the Keldysh formula,^® to 
describe the impact ionization transition rate. The Keldysh 
formula contains two parameters, P and which arc very 
difficult to determine from first principles. corresponds to 
the threshold energy for impact ionization while P is some 
measure of how strong the impact ionization transition rate is 
for caniers with energies greater than £n). A small value of 
P then corresponds to a soft threshold since the transition 
rate is relatively small near the threshold and increases only 
slowly thereafter. In contrast, a hard threshold is character- 
ized by a large value of P. Typically, both P and £, 1 , arc 
selected such that the calculated impact ionization rate is in 
agreement with experimental measurements. As such, more 
than one set of parameters can often reproduce the same 
ionization rate in a given material. Tang and Hess^ have 
shown that the calculated impact ionization rate in bulk sili- 


con, based on an ensemble Monte Carlo model including the 
Keldysh formula, agrees well with existing experimental 
data using either a soft or hard threshold model. 

Kane** made the first attempt to improve the impact ion- 
ization transition rate foimulation by directly evaluating the 
pair cross section using" pseudopotential-generated energy 
bands and wave functions Vor bulk silicon. In his formula- 
tion, Kane calculated the transition rate for different wave 
vector (k vector) states of comparable energy and averaged 
these rates to determine an energy dependent transition rate. 
More recently, Sano et and Bude and Hess^ have in- 
vestigated the threshold energy in silicon and GaAs using 
numerical techniques. Their work indicates that the threshold 
is highly k dependent and that the relative softness of the 
threshold originates from the k dependence of the impact 
ionization transition rate. 

New theories of interband impact ionization have been 
presented which treat the problem cither scmiclassically 
or quantum mechanically.^^'*^. The quantum-mechanical 
works have demonstrated that the threshold energy is always 
relatively soft owing to energy broadening effects arising 
from the electron-phonon interaction. In the quantum- 
mechanical formulations, a Lorentzian is substituted in place 
of the energy conserving delta function which leads to a 
broadening of the threshold state itself. Bude and 
co-workers*^ have demonstrated that through inclusion of 
quantum-mechanical broadening, a threshold state cannot 
rigorously be defined since the Lorentzian does not demand 
strict energy conservation over short time scales. As a result, 
a carrier can impact ionize even if its energy is less than the 
energy gap. Energy broadening acts then to inflate the impact 
ionization rate near threshold. 

SimUarly, even the semiclassical models predict that the 
ionization transition rate in bulk silicon is relatively low 
leading to a soft thresholdlike behavior. Owing to the fact 
that the ionization transition rate is calculated to be signifi- 
cantly lower than the competing phonon scattering rate, these 
models predict that the electrons drift to high energies, on 
average, before suffering an impact ionization event. As a 
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result, the uveragc energy of the electrons is expected to be 
quite high at high-electric-field strengths in contrast to what 
is predicted using a hard threshold model. Recent experimen- 
ttl measurements using x-ray photoemission sjKCtroscopy 
(XPS)’" indicate that the threshold for interband impact ion- 
ization in bulk silicon is exceedingly soft. These results for- 
ther indicate that the ionization rate may have multiple 
ranges of behavior which is consistent with theoretical 
predictions’’’^ that the second band dominates the impact 
ionization process in bulk silicon. ^ 

Aside from the obvious importonce of improving the un- 
derstanding of inteibahd impact ionization, knowledge of the 
nature of the threshold energy in different materials is useful 
in modeling advanced avalanche photodiode (APD) designs 
which utilize heterostructures and multiple quantum wells to 
enhaiKe the elefctron ionization rate.** In these devices the 
presence of the conduction-band-edge discontinuity formed 
at the heterostructure interface provides a kinetic energy 
boost to electrons injected into the narrow-band-gap layers 
from the wide-band-gap materials. The larger the con- 
duction-band-edge discontinuity is, the greater relative effect 
it has on the electron ionization rate. This depends though 
upon how soft the threshold is. If the threshold is relatively 
hard, the conduction-band-edge discontinuity can greatly in- 
fluence the ionization rate by providing sufficient energy to 
promote the electron from below threshold to above thresh- 
old. On the other hand, if the threshold is exceedingly soft, 
the presence of the conduction-band-edge discontinuity may 
not have as much of an effect since most of the carriers must 
attain very high energies before impact ionizing. As a result, 
the discontinuity may contribute 'relatively little to the over- 
all heating of the electron distribution. 

Since interband impact ionization is a threshold process, 
a carrier launched near the band edge, either initially or after 
suffering an ionization event, must drift for some time to 
acquire sufficient kinetic energy in order to impact ionize. 
The ionization probability then is necessarily equal to zero 
for some distance for both the generated carrier and the ini- 
tiating carrier immediately following an ionization event. 
The distance over which the ionization probability is as- 
sumed to be zero following launching from the band edge is 
often referred to as the dead space. The dead space can be of 
importance in multiquantum-well structures since if the well 
widths are much smaller than the dead space distance, then 
the ionization rate may be appreciably reduced, leading to 
little enhancement. Recently, the electron and hole ionization 
coefficients were measured in a series of GaAs/AlG^ i^l- 
tiple quantum wells with varying well and barrier widths. Jt 
was found that in the narrow well width structures (<100 A) 
there was virtually no enhancement of the experimentally 
measured electron to hole ionization rates ratio, while in 
wider well structures there was a significant increase in the 
ionization rates ratio. These results may be related to the 
dead space in the ionization process. 

In this article we present detailed calculations of the 
electron-initiated impact ionization rate in bulk Ga^ using a 
i jfTTiiHacticel formulation of the interband impact ionization 
transition rate. The transition rate is determined using Fer- 
mi’s golden rule from a two-body screened Coulomb inter- 


action assuming both energy and momentum conservation 
for the first two conduction bands of GaAs. These results are 
then incorporated into an ensemble Monte Carlo simulator to 
calculate the total impact ionization rate as a function of 
applied electric field in bulk GaAs. Aside from the ionization 
transition rale and the total impact ionization rate, the dead 
space for electron-initiated impact ionization is also calcu- 
lated and compared to results based on a relatively hard 
threshold model. The implications of an exaedingly soft 
threshold for impact ionization in GaAs arc discussed in re- 
lation to multiquantum-well, devices. The details of the 
model are summarized in Sec. II. The calculated results are 
presented in Sec. HI while conclusions arc drawn in Sec. IV. 

II. MODEL DESCRIPTION 

The of the model have been exhaustively re- 

viewed in a previous work.’* Nevertheless, for completeness, 
the salient features of our model arc briefly outlined here. 
The total impact ionization rate is determined using an en- 
semble Monte Carlo simulation which includes the full de- 
tails of the first two conduction bands of GaAs and the nu- 
merically generated semiclassical impact ionization transi- 
tion rate. The Monte Carlo simulator used in the c^culations 
is based on the original work of Shichijo and Hess, with the 
improved phonon scattering technique of Chang et al. and 
the improved band-structure interpolation technique of Fis- 
chetli and Uux.“ The material parameters used in the simu- 
lation arc taken from Refs. 27 and 28. The estimators used in 
the simulation are discussed at length in Ref. 28. 

The high-energy scattering rate is assumed to be domi- 
nated by deformation potential scattering. The deformation 
potential scattering rale is calculated using a full order time- 
dependent perturbation theory expansion assuming that tlw 
deformation potential is a constant, independent of energy. 
Although the assumption of a constant deformation potential 
may not be strictly valid,” it is standard and commonly used. 
The deformation potential scattering is determined by direct 
integration over the numerically generated density of states 
for both bands. At the points in k space at which both bands 
overlap, the deformation potential scattering depends, of 
course, on the density of stales within both bands. During the 
course of the Monte Carlo simulation, the electrons can 
transfer then from one band to the other through the action of 
the deformation potential scattering. 

Impact ionization is incorporated into the Monte Carlo 
model using the semiclassical, wave-vector (k)-dependent 
transition rate given by Ridley.^ The matrix element is as- 
sumed to be a screened Coulomb interaction between two 
electrons.’* The transition rate is evaluated numeriMlly by 
integrating oiler the foil Brillouin zone. Several million pos- 
sible final slates are sampled for each integration to ensure 
convergence and accuracy. The overlap integrals present in 
the matrix element are determined by diagonalizing a 15X15 
k-p Hamiltonian.” The values of the overlap integrals match 
closely to those presented by Burt et al.” The transition rate 
is determined for each mesh point within a finely spaced, 916 
point, k-space grid spanning the reduced zone * of the first 
Brillouin zone. Additional points are also used for interpola- 
tion purposes raising the total number of mesh points at 
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FIG. 1. CtlcuUtcd interband impact ionization transition rate for the first 
conduction band of bulk GaAs within the reduced zone of the first Brillouin 
zone. The ionization transition rate is plotted vs and ky bolding fixed at 
0.0. Those points which arc shown to have a transition rate equal to the 
minimum of the vertical axis scale have no threshold and hence have a zero 
transition rate. Owing to the logarithmic scale used in the figure, these 
points arc represented as having a minimum transition rate for convenience. 

I 

which the transition rate is calculated to 1419. The transition 
rate is determined for both conduction bands at each of the 
mesh points. The mesh employed for the impact ionization 
transition rate is the same as that for the energy interpolation 
used during the course of the Monte Carlo simulalion.^^ 

The instantaneous electron-impact ionization transition 
rate is determined both from its energy and k vector. When 
the simulated electron’s energy is less than the energy band 
gap, the electron cannot impact ionize and its flight is deter- 
mined in the usual manner.^^ When the simulated electron 
drifts to an energy greater than or equal to the energy band 
gap,' its ionization transition rate is determined using the k 
space mesh. The eight mesh points nearest to the components 
of the simulated electron’s k vector arc first determined. If 
the ionization transition rate is nonzero for all of the eight 
mesh points then the transition rate of the simulated electron 



0 0 


FIG. 2. Calculated interband impact ionization transition rate for the aecoiid 
conduction band of bulk GaAs within the reduced zone of the first Brillouin 
zone. The ionization transition nit is plotted vs k^ and ky holding k^ fixed at 
0.0. Those points which arc shown to have a transitioii rate equal to the 
minimum of the vertical axis scale have no threshold and hence have a zero 
transition rale. Owing to the logarithmic scale used in the figure, these 
points are represented as having a minimum transition rate for convenience. 



FIG. 3. Tblal calculated electron phonon scatierlng rate used in the loft 
threshold and hard threshold simulations. The high-energy scattering rate for 
the soft threshold case is matched at 0.5 cV with the phonon scallering rale 
calculated from Fermi's golden rule, while the rsle for Uie hard threshold 
case is matched at 0.8 eV. 


is determined by linearly interpolating between the values at 
each mesh point. If one or more of the mesh point values is 
zero, then the transition rate for the simulated electron is 
chosen equal to that of the nearest k point. The impact ion- 
ization transition rate determined in this manner is then 
added to the total phonon scattering rate corresponding to the 
energy of the simulated electron. As in standard Monte Carlo 
algorithms, a random number is thrown to determine if the 
carrier is scattered. Given that the carrier scatters, a second 
random number is generated to determine the particular 
event that occurs out of all the possible events at that energy 
and k vector. 

The calculated interband impact ionization transition rate 
for the first two conduction bands of bulk GaAs as a function 
of kj, and ky holding fixed at 0.0 within the reduced zone 
is presented in Figs. 1 and 2. As can be seen from compari- 
son of these two figures, the transition rate is significantly 
greater in the second conduction band than in 4hc first. The 
maximum transition rate in either band is only slightly 
greater than 10^^ s'’^ significanfiy less than the competing 
phonon scattering rale at these energies. The probability that 
an ionization event will occur at any given time, relative to a 
competing phonon scattering event, is then quite small. As a 
lesblt, an electron can drift for some time to very high energy 
before suffering an impact ionization event. 

III. CALCULATED RESULTS 

At high carrier energies the phonon scattering rate is 
assumed, as stated above, to .be dominated by deformation 
potential scattering. The scattering rate at these energies is 
determined by directly integrating over the final density 
states calculated from the numerically generated band struc- 
ture including collision broadening. A single deformation po- 
tential constant is selected in order to match the scattering 
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FIG. 4. Calculated impact ionization transition rates as a function of elec- 
Iron energy measured from the conduction-band minimum using the 
k-depcndcnl and Keldysh formula with two different sets of P and . The 
energy -dependent rate for the k-dependent rate is calculated by performing a 
further integration of the k-dependent rate. 


rate calculated from Fermi’s golden rule at a specific energy. 
Since the impact ionization coefficients are very sensitive to 
the total scattering rate, a small change in the total scattering 
rate will lead to a dramatic change in the impact ionization 
coefficients. Figure 3 shows the total electron phonon scat- 
tering rates with the deformation potential scattering rate 
matched to the scattering rate calculated from Fermi’s golden 
rule at 0.5 and 0.8 eV, respectively. As is discussed later, the 
lower scattering rate, that matched at 0.5 eV, results in a 
softer threshold for the impact ionization process, while the 
higher scattering rate, that matched at 0.8 cV, results in a 
harder threshold. It should be noted that there is some arbi- 



HG. 5. Calculated and experimental results for the electron-impact ioniza- 
tion coeffidcnta in bulk GaAs plotted as a function of inverse electric field, 
using the k-dependeni model, the Keldysh formula model with F=0.02 and 
£rt,"1.85 eV, and the Keldysh formula model with F«0.2 and E^* 1.42 cV. 
The expcrimenlal results arc taken from Bulman ct qI (Ref. 33). 
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FIG. 6. (a) Histogram ahowing the percentage of electron ionization events 
originating from the first and second conduction bands as a function of 
applied electric field based on the k-dependent model, (b) Histogram show- 
ing the percentage of electron ionization events originating from the first and 
second conduction bands as a function of applied electric field based on the 
Keldysh formula model with P*0.02, eV. (c) Histogram showing 

the percentage of electron ionization events originating from the first and 
second conduction bands as a function of applied electric field based on the 
Keldysh formula model with F^O.2. Ed, » 1.42 eV. 
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trariness in the total scattering rate since the deformation 
potentials are not fully known as a function of energy. Topi- 
cally, the scattering rate is adjusted such that the calculated 
canier drift velocities agree with the experimental measure- 
ments; however, the drift velocity is less sensitive to the 
scattering rate than the impact ionization coefficients. In 
these calculations, both sets of scattering rates show good 
agreement with the experimental saturation carrier drift ve- 
locity within the error tolerance of the Monte Carlo tech- 
nique. Therefore, either total scattering rate is plausible. 

An energy-dependent impact ionization transition rale 
can also be derived from the k-dependent rate by performing 
an additional integration over the initiating electron k-vector 
space, normalized by the density of states. The result is plot- 
ted as a function of initiating electron energy in Fig. 4, along 
with the rates using the Keldysh formula with two sets of P 
and £=0.02 and £t|j=l,85 cV, as well as £=0.2 and 
£ijj=1.42 cV. Multiple sets of £ and £,h may be found by 
fitting the calculated impact ionization coefficients to the ex- 
periment. The first set of values of £=0.02 and £,^=1.85 eV 
is chosen so as to match the Keldysh rate as close as possible 
to the energy-dependent rate derived from the k-dependent 
rate. Tlie value of 0.02 for £ is in the range of values which 
typically corresponds to a soft threshold condition. Another 
set of parameters for the Keldysh formula, £=0.2 and 
£(j,=1.42 eV, is chosen to compare the difference betwden 
the hard and soft threshold models. 

The impact ionization coefficients calculated using the 
k-dependent impact ionization transition rate, the Keldysh 
formula with £=0.02 and £ih=1.85 eV, and the Keldysh 
formula with £=0.2 and E^^=\A2 eV along with the experi- 
mental measurements of Bulman et al?^ arc plotted in Fig. 5 
at an electric-field direction of [100]. The lower total scatter- 
ing rate as shown in Fig. 3 is used in the calculations of both 
the k-dependent model and the Keldysh formula with 
£=0,02 and £jj,=1.85 eV. In the case of the Keldysh for- 
mula with £=0.2 and £j|,=l,42 cV the higher total scatter- 
ing rate must be used due to the higher impact ionization 
transition rate and the harder threshold present in this case in 
order to agree with the experimental results. As can be seen 
from Fig. 5, the k-dependent calculations match the experi- 
ment fairly well throughout the range of applied electric-field 
strengths plotted here. This is excellent considering there is 
no adjustable parameter used in the model for the impact 
ionization transition rate. The only ambiguity in this case 
appears in the scattering rate. The calculated results using the 
Keldysh formula for the impact ionization transition rate, 
within the otherwise identical Monte Carlo model, are in- 
cluded for comparison. With two adjustable parameters, £ 
and £(h it is possible to have better agreement with the ex- 
periments than the k-dependent model, as shown here for the 
case of £=0.02 and £,|,=1.85 cV. The results calculated us- 
ing the k-dependent model are closer to the results calculated 
by the Keldysh formula with a softer threshold as expected, 
since the k-dependent model also has effectively a very soft 
threshold. 

The effect of the second conduction band on the impact 
ionization coefficients in bulk GaAs for all three impact ion- 
ization transition rate models has also been investigated. Fig- 


ure 6(a) shows the percentage of electron ionization events 
originating from the first and second conduction bands 
within the k-dependent model. Nearly 100% of the impact 
ionization events originate from the second conduction band 
due to the fact that the impact ionization transition rate is 
much higher in the second conduction band. The percentage 
of electron ionization events originating from the first and 
second conduction bands for both Keldysh formula models 
are plotted in Figs. 6(b) and 6(c). In contrast to the 
k-dependent model, more ionization events originated from 
the first conduction band in both cases. This can be under- 
stood since in the Keldysh formula, the impact ionization 
transition rate is isotropic and no distinction is made between 
states in the first and second conduction bands. Therefore, 
the electrons have the same ionization probability at a given 
energy from the first conduction band as from the second 
conduction band in the Keldysh formula model, while in the 
k-dependent model the ionization probability is typically 
higher in the second conduction band than the first conduc- 
tion band even at the same energy. Hence, it is not surprising 
that the Keldysh formula model predicts that more events 
will occur from the first conduction band than the 
k-dependent model does. As the threshold becomes lower 
and harder, more ionization events occur in the first conduc- 
tion band. Therefore, the Keldysh formula model will predict 
very different physical results if a different set of £ and 
is used. 

It is interesting to examine the impact ionization coeffi- 
cients along different applied electric-field directions. The 
directional dependence of the impact ionization transition 
rate is incorporated naturally in the k-dependent formulation. 
Calculations of the electron-impact ionization coefficient at 
different applied electric field directions of [100], [110], and 
[111] show no field directional dependence of the impact 
ionization coefficients even with the k-dependent impact ion- 
ization transition rate. These results agree with previous cal- 
culations by Brennan and Hess^ using the Keldysh formula 
and the calculations reporled here using both Keldysh for- 
mula models. As shown in Fig. 6(a), almost all the impact 
ionization events originate from the second conduction band 
using the k-dependent formulation. Since the impact ioniza- 
tion transition rate in the second conduction band, shown in 
Fig. 2, is higher and more uniform than the rate in the first 
conduction band, little electric-field directional dependence 
of the electron-impact ionization coefficients is expected. 

The electron number density functions, defined here as 
the product of the electron density of states function and the 
electron distribution function, arc plotted in Figs. 7(a)-7(c) 
for all three models at an applied electric field of 250 kV/cm. 
The number density functions at the applied electric field of 
500 kV/cm for all three models are presented in Figures 
8(a)~8(c). Both the k-dependent model and the Keldysh for- 
mula model with the softer threshold predict similar electron 
distributions. The second peak in the distributions at an elec- 
tric field of 500 kV/cm for the k-dependent model and the 
Keldysh formula with softer threshold indicate many elec- 
trons occupy the second conduction band at high applied 
electric fields. In the Keldysh formula model with the harder 
threshold, the distribution is much cooler. In this case, fewer 
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FIG. 7. (a) The dccuon number function, defined as the product of the 
electron dcnsity of-slalcs function and the electron distribution function, at 
the applied electric field of 250 kV/cm calculated using the k-dependent 
model. Notice the significant number of electrons within the high-energy tail 
of the distribution. In the k-dependent model owing to the ultrasoft thresh- 
old present, the impact ionization process docs not truncate the high-energy 
tail of the distribution, (b) The electron number function, defined as the 
product of the electron dcnsiiy-of-stalcs function and the electron distribu- 
tion function, at the applied electric field of 250 kV/cm calculated using the 
Keldysh formula model with F-0,02, £<,-1.85 eV. As in the k-dcpcndcnl 
model, the high-energy tail of the electron distribution is not truncated due 
to the soft threshold in the ionization process, (c) The electron num^r 
function, defined as the product of the electron density-of-states function 
and the elecUon distribution function, at the applied electric field of 250 
kV/cm calculated using the Keldysh formula m^cl with r-0.2, 
eV. Notice that few electrons exist within the high-energy tail of the distri- 
t button in this case. The presence of a hard threshold acts to sharply truncate 
the high-energy tail at relatively low carrier energy. 



no. 8. (a) The electron number function, defined u the product of the 
electron density-of-states function and the electron distribution function, at 
the applied electric field of 500 kV/cm calculated using the k-dependent 
model. The second peak in the distribution is due to carriers within the 
second conduction band. Notice that there is a significant population of 
electrons at very high energy due to the presence of a very soft threshold for 
impact ionization, (b) The electron number function, defined as the product 
of the electron density-of-states function and the electron distribution func- 
tion. at the applied electric field of 500 kV/cm calculated using the Keldysh 
formula model with P-0.02, £*-1.85 eV. As in the k-dependent case, the 
second peak in the distribution is due to carriers within the second conduc- 
tion band. Many electrons survive to high energy in this case owing to the 
very soft threshold for impact ionization, (c) The electron numlxr function, 
defined as the product of the electron dcnsity-of-slalcs function and the 
electron distribution function, at the applied elcctnc field of 500 kV/cm 
calculated using the Keldysh formula model with P-0.2, £*-1.42 eV. 
Owing to the relatively hard threshold used in this case, the high-energy tail 
of the distribution is effectively truncated. 
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G 9 (•) Scalier plot of electron energy and the location at which impact 
nization occurs at the applied electric field of 333 kV/cm based on the 
dependent model. The dead space can be estimated by the distance mea- 
ired from lero to the first ionization event. In this case, this distance is 
90.0 nm. Notice that the average energy at which an im^ct ionization 
rent occurs is quite high, always greater than 2 cV and in inany cases 
loscr lo 3-4 cV. (b) Scatter plot of electron energy tnd the location at 
hich impact ionization occurs at the applied electric field of 333 kV/OT 
ased on the Keldysh formula model with 0.02, £#,-!.« eV. ^ in the 
.dependent case, the average energy at which an impact ionization event 
ccurs is relatively high and the dead space is relatively long, (c) Scatter 
lot of electron energy and the location at which impact km^lion ^rs a 
le applied electric field of 333 kV/cm based on the Keldysh formula mt^el 
rith £=0.2, £*-1.42 eV. Notice that the average energy at which an im- 
act ionization event occurs is relatively low. Many ionizations occur at 
nergies less than 2 eV and few carriers survive to energies much greater 
lian 2 eV. Notice also that the dead space is predicted to be much less than 
n the k^dependent model. In this case, the dead space is -50.0 nm almost 
lalf of that predicted for the k-^Jependenl model. 


electrons survive to high energies. This is, of course, due to 
the fact that in the hard threshold model, once a carrier at- 
tains threshold it is highly likely to suffer an impact ioniza- 
tion event. As a result, few carriers drift to energies signifi- 
cantly higher than the threshold. 

The electron energies and positions at which impact ion- 
izations occur are plotted in Figs. 9(a)-9(c) for all three 
models at the applied electric field of 333 kV/cm. Each dot 
in the figures represents an impact ionization event that oc- 
cured at energy £ and distance D from launch. The electrons 
arc all launched initially from zero energy. From inspection 
of the figures, the electron-impact ionization dead space can 
be estimated at each specific applied electric field. Notice 
that the dead space predicted by the harder threshold inodel 
is the shortest among the three models. This is because m the 
harder threshold model electrons can impact ionize at much 
lower energies [Fig. 9(c)] due to the lower threshold and the 
higher impact ionization probability once the carrier reaches 
the threshold. The distance traveled for the electrons to reach 
the ionization energy thus is much shorter. The k-dependent 
model has the longest dead space. In the k-dependent model 
most impact ionizations occur in the second conduction 
band, while in the Keldysh formula model electrons can 
more easily impact ionize within the first conduction band. 
Since electrons can only Uansfer to the second conduction 
band through scattering events in our simulation model, it 
takes a longer time and, thus, a longer distance, for rhe elec- 
trons to impact ionize. Due to the use of a higher and softer 
threshold, electrons need reach higher energy to impact ion- 
ize. This results in longer dead spaces than that reported in 
Ref. 35. 

IV. CONCLUSIONS 

In this article we have presented an analysis of the inter- 
band impact ionization rate in bulk GaAs using a semiclas- 
sical, k-dependent, numerical formulation of the impact ion- 
ization transition rate. It is found that the transition rate is 
substantially higher within the second conduction band than 
the first conduction band. As a result, the overwhelming ma- 
jority of electrons impact ionize from states within the sec- 
ond conduction band than from the first conduction band, in 
contrast to what is predicted by a hard threshold, Keldysh 
formula model. In addition, the numerical, k-dependent for- 
mulation of the transition rate predicts that the electrons sur- 
vive to very high energy, on average, before suffering an 
impact ionization event. Again, this is vastly different from 
that predicted by the hard threshold, Keldysh formula model 
wherein most of the electrons ionize at much lower energies. 

Interestingly, both the hard threshold model, and the soft 
threshold models, i.e., either a soft threshold, Keldysh for- 
mula model or the numerical, k-dependent model, can be 
fashioned, by adjusting the total phonon scattering rate, to 
match the experimental results exceedingly well. Though the 
hard and soft threshold models require different total phonon 
scattering rates in their description, the two different rates 
examined here both lead to acceptable predictions of the 
steady-state electron drift velocities. As such, it is difficult 
under present circumstances to definitively choose which 
scattering rate is the most appropriate. A more detailed 
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analysis of the high-energy scattering rate, taking into ac- 
count the energy dependence of the deformation potentials, 
needs ultimately to be done. Nevertheless, it is expected that 
the k-dependent transition rate formulation is far more accu- 
rate than the Keldysh formula. In order that the k-dependent 
transition rate calculations agree with the experimental mea- 
surements, the phonon scattering rate would have to be lower 
than that usually quoted for a hard threshold model. 

If the threshold for interband impact ionization is indeed 
very soft in GaAs, as indicated by the present k-dependent 
transition rate calculations, the dead space for impact ioniza- 
tion is then significantly larger than previously estimated. A 
longer dead space may lead to different predictions of the 
electron-impact ionization enhancement in GaAs/AlGaAs 
multiquantum-well (MQW) structures. As a result, the opti- 
mal well and barrier widths of these devices may be quite 
different from those predicted using a harder threshold 
model. We are currently investigating how an ultrasoft 
threshold, as predicted by the k-dependent transition rate cal- 
culations presented here, would effect the predicted electron 
and hole ionization rales in MQW structures. These results 
will be presented elsewhere. 
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TYieoretical Investigation of wave-vector-dependent analytical 
and numerical formulations of the Interband Impact-lonizatlon 
transition rate for electrons In bulk silicon and GaAs 
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(Received 21 March 1994; accepted for publication 8 June 1994) 

Tbe electron interband impact-ionization rate for both silicon and gallium arsenide is cal^lated 
using an ensemble Monte Carlo simulation with the expressed purpose of comparing different 
formulations of the interband ionization transition rate. Specifically, three different UMtmenIs of the 
transition rate are examined: the traditional Keldysh formula, a new k-de^ndcnt 
formulation first derived by W. Quade, E. Scholl, and M. Rudan [Solid St«‘e Ele^on. 36 1493 
(1993)], and a more exact, numerical method of Y. Wang and K- F. Brennan [J. Appl. Phys. 75, 313 
(1994)]. Although the completely numerical formulation contains no adjustable parameters and as 
such provides a very reliable result, it is highly computationally intensive, ^ternativel^ the Keldysh 
formula, although inherently simple and compuUtionally efficient, fails to include the k dependent 
as well as the details of the energy band structure. The k-dependent analytical formulation of Quade 
and co-workers overcomes the limitations of both of these models but at the ex^nse of some new 
parameterization. It is found that the k-dependent analytical method of Quade and ^workers 
produces very similar results to those obtained with the completely nument^l model for some 
quantities. Specifically, both models predict that the effective threshold for impact ionization in 
GaAs and silicon is quite soft, that the majority of ionization events onginate f^ the second 
conduction band in both materials, and that the transition rate is k dependent. Therefore, it 
concluded that the k-dependent analytical model can qualitatively reproduce results similar to Aose 
obtained with the numerical model yet with far greater computational efficiency. Nevertheless, there 
exist some important drawbacks to the k-dependent analytical m^el of 
These arc that it does not accurately reproduce the quantum yield data for bidk 
determination of a new parameter, related physically to the overlap integrals of the Bloch state 
which can only be adjusted by comparison to experiment, and fails to account for any wave-vector 
dependence of the overlap integrals. As such the transition rate may be overestimated at those poinU 
for which “near vertical,” small change in k, transitions occur. 


I. INTRODUCTION 

The theoretical study of inlerband impact ionization has 
been greatly aided by the advent of numerical methods. C^- 
ing to the complexity of the encigy band structure at high 
energies, at which impact ionization is initiated, w well as 
the complicated carrier-phonon scattering mechanisms, sim- 
plified analytical formulations, such as the early theories of 
Shockley,' Wolff,^ and Baraff^ have limited validity. Nu- 
merical models for calculating the ionization rate were then 
advanced based on the Monte Carlo method; however, 
these approaches utilized parabolic or nonparabolic analyti- 
cal energy bands which are of questionable validity at high 
carrier energies. The full details of the energy band structure 
were first accounted for in the Monte Carlo model of 
Shichijo and Hess.* In their model the dynamics of the elec- 
trons in bulk GaAs were simulated within the first conduc- 
tion band, calculated based on an empirical pseudopotential 
model. Later, Fischelti and Laux’ developed a more ad- 
vanced Monte Carlo simulator for studying impact ionization 
which improved the numerical accuracy as well as incorpo- 
rated transport within higher conduction bands. 

In order to calculate the overall impact ionization rate 
using the Monte Carlo technique, it is necessary to formulate 


an expression for the impact ionization transition rate. This 
formulation is then incorporated into the Monte Carlo simu- 
lator and is treated as an additional scattering mechanism. An 
impact-ionization event is chosen stochastically in the usual 
way* through the use of a random number based on the rela- 
tive magnitude of the ionization transition rate compared to 
competing phonon scattering events. Although there has 
been recent work on formulating the ionization transition 
rate using a higher-order quantum-mechanical perturbation 
theory expansion,*”" the impact-ionization transition rate is 
more typically determined from use of Fermi’s golden rule. 
Fermi’s golden rule can be expressed as*^ 

Wn- J ^ \M\^S{Ef~E,)S{mdSf, (D 

where dSf represents integration over all the final states, E, 
and Ef are the energies of the initial and final states, res|»c- 
tively, Ak is the minnentum change during the interaction. 
The matrix element M in Eq. (1) is formed using toth the 
direct and exchange M, terms following Ridley as 

( 2 ) 

where 
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/(k,.k;)/(kj.ki) 

€V |k;-k,|^+X^ ’ 

(3) 

e* /(k,.ki)/(k2.k:) 

€V jkj-kil^+X^ ■ 

In Eq. (3), kj ,k{ arc the wave vectors of the incident elec- 
tron before and after the interaction. k|, kj represent the 
states of the secondary electron-hole pair after the collision, 
X is the static screening factor, and e is the electronic charge; 
V is the crystal volume, and e is the dielectric function. 
Quantities /(k,,k;),/(kj,ki),/(k,,ki),/(k 2 k;) arc overlap 
integrals.’^ Following Ridley,** the term |A/^— A/,|* is ne- 
glected and the direct and the exchange terms are assumed to 
be equal. 

In roost previous Monte Carlo simulations (e.g., Refs. 6 
and 13) the impact-ionization transition rate is calculated us- 
ing the Keldysh formula,*^ which can be obtained from Eq. 
(1) by adc^ting several simplifying assumptions: These are 
that the transition is calculated assuming a direct semicon- 
ductor with parabolic bands for all of the carriers, the overlap 
integrals are considered to be constant and the denominator 
in the matrix element in Eq. (3) is also taken as a constant 
with thresholds values of kj ,k| . This yields the very well- 
known quadratic dependence of the transition rate Wy on the 
energy E of the initiating particle, usually expressed as 

where E^, is the threshold energy and W(£rt,) is the total 
phonon scattering rate at threshold. £^, and the prefactor p 
cannot generally be determined from first principles; instead, 
they are determined through comparison of the calculations 
to experimental data. Aside from this obvious limitation of 
the Keldysh formula, the neglect of an accurate accounting 
of the high-energy region of the band structure also results in 
the failure to properly assess the k dependence of the ioniza- 
tion transition rate itself. Recent theoretical studies have in- 
dicated that there is a significant deviation in the ionization 
uansitioD rate with respect to the initiating carrier’s k vector 
in different roaterials.*^’*‘ 

Several alternative approaches to the Keldysh formula 
for calculating the interband impact ionization rate have been 
recently presented.*^'** In each of these methods, the inter- 
band inqMct-ionization transition rate is evaluated numeri- 
cally and can be then incorporated into an ensemble Monte 
Carlo calculation. Although these approaches presently offer 
a far more accurate treatment of impact ionization than the 
Keldysh formula, they require extensive numerical computa- 
tion. In the approach adopted by Wang and Brennan*^'*^ the 
transition rate is determined directly from Fermi’s golden 
rule by numerically integrating Eq. (1) at different k points 
within the first Brillouin zone. Their analysis includes the 
direct calculation of the overlap integrate from the numeri- 
cally generated wave functions using a k*p method. The 
transition rate is evaluated by integrating over several mil- 
lion final states directly evaluating the overlap integrate in 
each case. This is particularly computationally intensive 
since the overlap integrals must be evaluated for each final 


state. Sano and Yoshii'* have used a less computationally 
intensive scheme. In their method, the three-dimensional k 
space incorpomting the reduced zone is partitioned into 
«m«il cubes. The matrix elements within each cube are as- 
sumed to be equal alleviating the need to evaluate the over- 
lap integrate for each possible state. In this way, the number 
of computations can be greatly reduced at the expense of 
increased memory r^uirements. However, we have found 
that the overlap integrals vary in a nonpredictive manner 
from point to point in our k*p calculations and cannot gen- 
erally be assumed to be the same even for states very close in 
k space. For this reason and memory storage limitations, in 
the model of Wang and Brennan”’** the overlap integrate are 
evaluated for every final state. Nevertheless, all of the nu- 
merical approaches indicate that the transition rate is depen- 
dent on ^ initiating electron wave vector. 

Recently, another analytical treatment of expression (1), 
which also includes the k-vector dependence of the impact- 
ionization transition rate, was presented by Quade, Scholl, 
and Rudan.** They calculated the impact-ionization transi- 
tion rate for the general case of an arbitrarily shaped energy 
band for the inqract-ionizing conduction electron, and three 
anisotropic parabolic bands for the final states with their ex- 
trema located at different points in the k-vector space. The 
assumption of the parabolic energy bands for the final states 
is acceptable since the final states are generally at low energy 
near the band minimum where the parabolic approximation 
is typically good. Other approximations used in the deriva- 
tion include nondegeneracy and a constant value for the 
overlap integrals. Although the formulation developed by 
Quade and co-workers** contaitu some parameterization due 
to the tact that the overlap integrals are not evaluated but are 
assumed to be constant parameters, it offers a far more com- 
putationally efficient means of evaluating the ionization tran- 
’ sition rate with the inclusion of the k dependence than the 
numerical method of Wang and Brennan.**'** Subsequently, 
it is of interest to examine bow the k-dependent, analytical 
model of Quade and co-workers** compares to the numerical 
model of Wang and Brennan.**’** 

In this article, we incorporate the formulation of the 
impact-ionization transition rate of Quade and co-workers** 
within an ensemble Monte Carlo simulator to determine the 
electron impact-wnizatkm rate in bulk silicon and GaAs. The 
impact-ionization rate calculated in this way is then com- 
pared to experimental data as well as calculations made in an 
otherwise identical Monte Carlo simulator using the Keldysh 
formula and the direct, k-dependent transition rate. The na- 
ture of the ionization rate is further probed using the Quade 
and co-workers’ model to determine if it yields similar physi- 
cal results to that of the ilirect, numerical model. Specifically, 
the transition rate udthin eadi band is calculated and the 
percentage of ionization events originating from each con- 
duction band is determined. The number density distributions 
and the quantum yield as a function of energy are determined 
and compared between the ilifferent theoretical models as 
well. In Sec. II the details of the theoretical models are re- 
viewed. The calculated results are presented in Sec. Ill and 
conclusioiu are drawn in Sec. IV. 
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II. MODEL DESCRIPTION 


The total impact-ionization rate as a function of inverse 
electric field is calculated using an ensemble Monte Carlo 
simulation, based on the original work of Shichijo and Hess.^ 
The details of the simulator have been extensively described 
elscwhcre;*^’*^ therefore, only the main features arc outlined 
here. The simulation includes the full details of the first and 
second conduction bands for silicon or gallium arsenide. The 
scattering rate in either case is determined using the im- 
proved phonon scattering technique of Chang et al?^ and the 
improved band-structure interpolation technique of Fischetti 
and Laux.^ The total phonon scattering rate is determined in 
the following way. Within the low-energy region, the scatter- 
ing rate is calculated from Fermi’s golden rule for all the 
relevant mechanisms present.^^ Within the high-energy 
range, deformation potential scattering is assumed to be the 
dominant scattering mechanism. The scattering rate b this 
region is obtained by integrating over the final density of 
states, calculated from the numerically generated band struc- 
tures, including collision broadening of the final state. The 
deformation potential is assumed constant and is selected to 
match the scattering rate calculated from Fermi’s golden rule 
at a low specific energy. Interband electronic transitions arc 
enabled via the action of the deformation potential scattering 
as described in Ref. 15. The actual scattering rates used in 
these calculations have been published previously in Ref. 15 / 
for silicon and Ref. 16 for GaAs. 

The impact-ionization transition rate is formulated in 
three different ways and is then incorporated into the Monte 
Carlo simulator described above. The different transition rate 
formulations used are a direct, numerical method previously 
described in detail in Refs. 15 and 16, the Keldysh formula 
described in detail in Refs. 6 and 7, and a new k-dependent 
analytical formulation derived by Quade and co-workers^ 
which is summarized below. 

The k-dependent analytical formulation of the impact- 
ionization transition rate, hereafter referred to as the Quade 
formula, predicts the ionization transition rate for an initiat- 
ing electron of wave vector kj , colliding with another elec- 
tron in the valence band, to bc^ 

where 


ifl /£, 


(6) 

1 1 >/m I 


(7) 

T“4^nf.,(a,+A)''M 

t 


“ /V fl 

«i+pi 


,(8) 


g(k,)=^ |£o+2 


( 9 ) 

( 10 ) 


and 

In the above set of equations, a, fi, and y represent a para- 
bolic vdence band and two possibly different parabolic con- 
duction bands for the final states of the hole and electrons 
after the impact ionization event, respectively. The aniso- 
tropic bands are parametrized by ratios of effective masses 

where i labels 

the Cartesian coordinate, and the extremal points of each 
band are positioned at Z., Z^, Z,; is an isotropic effec- 
tive mass of the conduction band for the initiating electron. F 
represents the product of the squared overlap integrals and is 
taken as constant for each, initiating band but may vary be- 
tween bands, /t stands for the geometrical average of the 
quantities /i|, fi 2 , defined by Eq. (8). The arbitrarily 
shaped function £(k|) represents the energy of the initiating 
electron, Eq represents the gap plus offsets of the extremal 
points of the bands a, fi, y. Iq“k|-L, where L is a recip- 
rocal lattice vector and 6 is the Heaviside function. Equation 
(5) is derived assuming a static dielectric constant. Although 
recent work*’ has shown that a full wave-vector, frequency- 
dependent dielectric function can modify the ionization tran- 
sition rate, its incorporation in the Quade formula would not 
lead to an analytic^ expression. For this reason, a constant 
dielectric function is used within the Quade formulation. To 
facilitate comparison of the Quade model to other models, a 
constant dielectric function is used within the numerical and 
Keldysh formula models as well. 

The total impact-ionization scattering rale for the initiat- 
ing particle with wave vector k| is calculated taking into 
account all of the possible final states within the full Bril- 
louin zone. The total impact-ionization transition rate can 
then be expressed as 


»Fu(k,)“2 2 2 Wg(k,)..p.,, (12) 

^ fi y 

where the summation goes over all the valence bands a and 
all possible valleys y in the first and second conduction 
bands. The summations in Eq. (12) arc taken so as to include 
both direct and exchange terms. 

The material parameters used in the calculations have 
been taken from Refs. 25-28 for silicon and Rcfe. 29-32 for 
gallium arsenide. For the case of the second conduction band 
in both materials, the longitudinal m| and transverse ef- 
fective masses for the X valleys have been obtained from the 
numerically calculated band structures determined by a k*p 
calculation. The longitudinal and transverse masses used for 
silicon and gallium arsenide in the second conduction band 
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TABLE 1. Eneifie* of tymmeuy points (relative to the top of the valence 
band, in eV) and values of relative effective masses used in the calralations 
of the impact-iooisatioo transition rate in silicon and gallium arsenide. 

Silicon Gallium anenide 


E(Xa) 

E(Lci) 

E{Xc„) 

E(Tso) 

m^Xci) 

m^Xa) 

m^Xc) 

oiAXcii) 

"(r Cl) 
"(rnn) 
■•(rtn) 
"(rjo) 


3.4? 

I.l? 

1.9? 

1.14' 

-0.044‘ 

Ol? 

0.9? 

0 . 21 * 

0.29* 

0.284* 

0.537^ 

0.153' 

0.234' 


1.423" 

1B98" 

1.706" 

1.973' 

-0.340" 

0.27* 

1.98' 

0.23* 

0.40* 

0.18' 

0.06f 

0.45' 

0.082' 

0.170' 


‘Reference 32. 

^Reference 30. 

*\%lues obtiined from the oumerically calculated band structure. 
^Reference 26. 

‘Reference 25. 

'Reference 31. 


arc reported in Table I. In the case of L valleys, F valleys, 
and all the valence bands, the isotropic effective masses have 
been used and their values are also listed in T^ble I, as well 
as the energies of symmetry points relative to the top of the 
valence band. 

The inteiband impact-ionization transition rale calcu- 
lated using the approach based on Quade and co-workers 
results assuming that the overlap integrals are simply unity is 
presented in Figs. 1(a) and 1(b) for the first and second con- 
duction bands in silicon, and in Figs. 2(a) and 2(b) for the 
first and second conduction bands in gallium arsenide, re- 
spectively. Comparison of Figs, 1 and 2 to similar figures 
obtained by Wang and Brcnnan*^-’^ using the completely nu- 
merical technique in which the overlap integrals are evalu- 
ated directly shows good qualitative agreement between the 
two approaches. For bulk silicon the primary difference, 
aside from the overall magnitude, is that the number of 
points for which a nonzero transition rate has been calculated 
is higher using the completely numerical approach. A similar 
effect can be seen in the case of the first and second conduc- 
tion bands in gallium arsenide. Again notice that there is a 
significant difference in the number of nonzero transition rate 
points between the two models (data calculated using the 
numerical approach can be found in Ref. 16) for the first 
conduction band; however, the rate from these points is rela- 
tively small compared to the competing scattering rates. 

Further examination of Figs. 1 and 2 reveals that the 
transition is greater for electrons originating from the second 
conduction band than from the first conduction band, which 
corresponds to the results based on the numerical 
calculations.*^'^^ This is particularly true for G^. As is 
discussed below, due to the fact that the transition ^ 
greatest for the second conduction band, most ioniMtion 
events occur from electrons which originate from within the 
second conduction band. 



FIG. 1 , CalcuUled inleiband impaci-ionteaiion tniuttion rile for the firsi 
tnd second conduction bind in bulk silicon is 1 function of snd ky (k, is 
fixed It 0.0) within the reduced xonc of the first Brillouin zone obuined by 
the inslyticil formuUiiwi based on the Quide model. InitiiUng electron is 
in (1) the first conduction bind and (b) the second conduction band. 


III. RESULTS 

The impact-ionization rate as a function of inverse elec- 
tric field, calculated using the ensemble Monte Carlo simu- 
lation, is presented in Figs. 3 and 4 for silicon and gallium 
arsenide, respectively. In each case, the rate is calculated 
using the Keldysh formula, the k-dependent numerical for- 
mulation, and the Quade formula for the impact-ionization 
transition rate. The experimentally obtained dau are also 
shown in both figures. The experimental results of van Over- 
atraeten and DeMan,^ Grant, ^ and Woods, Johnson, and 
Lambert^ for electron-initiated impact ionization in silicon 
are plotted in Fig. 3. These data comprise a represenUtive set 
of data which cover the full range of experimental measure- 
ments for electron-impact ionization in silicon. The experi- 
mental measurements of Bulman et are plotted for gal- 
lium arsenide in Fig. 4. 
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FIG 2 C»lcuUled inletbind impacl-ioniz»lion Uinsition rtle tor Ihe first 

and second conduction band in bulk GtAs as a function of *, and *, (*. is 
fixed at 0.0) within the reduced wne of the first Brillouin tone obtained by 
the analytical formulation based on the Quade model. Initiating el^ron is 
in (a) the first conduction band and (b) the second conduction band. 


In the case of silicon, good agreement with the daU of 
Woods and co-workers is obtained for the k-dependent to- 
merical formulation, as was also reported previously. No 
adjustable parameters for the impact-ionization transition 
rate arc used in this case. Very similar results can be obtained 
yeine the Quadc formulation by choosing the values for Uie 
parameter F in Eq. (7) equal to 0.24 for both the first and 
second conduction bands. F corresponds to an avera^ 
value of the squared overlap integrals. The squared overlap 
integrals always have value between 0 wd 1. In the analyU- 
cal model of Quadc and co-workers and that presented 
here, F is treated as a parameter which is chosen so as to 
yield good agreement with experiment. Best agreement with 
the low range of silicon experimental data occurs for 
F*=0.24 given the phonon scattering rate adopted here. The 
calculated silicon impact-ionization rate obuined using the 



FIG. 3. OilcuUled awl experimenul rewltt for the election hnpirt- 

ionizalion eoefBcieoU in bulk xllicon ploMed u i fbnetion 

field. Retolls of rhe eriaiUliom of the tonlMtion nie indoding the Keldyih 

fomula, k-dep e ndent mmierical, and Quade model are preien^d. Expen - 

mental ftaoha of Refc. 33 , 34, and 35 are •**o J*own 

the values of the averag^ squared overtop integiab used In the Quade 

• • « a •" 


Keldysh formula with p*0.02 and £^,*=1.13 eV is also 
shown in Fig. 3. A value of p of about 0.02 is typically 
ooiisidered to correspond to a soft threshold. In contrMt to 
both k-dependent models, belter agreement with the higher 
range of experimental measurements can be seen for this 
particular choice of threshold and p; however, a different 
choice of the parameter F within the Quade formulation 
would also provide a better fit to the higher range of experi- 
mental data as well. 

In gallium arsenide the dependence of the impact- 
ionization rate on the inverse electric field obtained using the 
numerical, the Quadc models with F=0.05 and 0.005 for 
the first and second conduction band, respectively, and the 
Keldysh formula with p=0.02, £o,= l-85 eV for the 
impact-ionization transition rate, shows in all three cases 
fairly good agreement with the experimentally measured data 



FIO. 4. Ctolcutoled and cxperimcntol leaulu fot the etectrw tapact- 
iontoaiion coefficknto in bulk griUum anenide plotted •$ a function of 
Imytne electric field, using the KeWyah formula, ’ 

and Quade model. Experimenul leaulta aiiown are those of Bulitun « al. 
(Ref. 36). FI and F2 represent ihe values of the averaged squai^ overtop 
integrals used in the Quade model tor bands 1 and 2, respective y. 
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FIG. 5. Hislogram showing the percemage of ionization events originating 
from the first and aecood conduction bands in silicon for (a) the k-dcpendent 
numerical formulation and (b) the Quade formulation. 


(Fig. 4). It should be noted that in the case of the approach 
based on the Quade formula, the selection of the fitting factor 
F for the first conduction band essentially does not play any 
role due to the negligible impact-ionization transition rate of 
the electrons in the first conduction band with respect to the 
second conduction band. 

The percentage of ionization events originating fiom the 
first and the second conduction band in bulk silicon as a 
function of applied electric field is plotted in Fig. S(a) for the 
numerical model and in Fig. 5(b) for the model based on the 
Quade formula. At the applied electric fields investigated 
here the impact-ionization events originate predominantly in 
the second band in both cases. Therefore, on the basis of 
these reports, it appears that most electron-initiated ioniza- 
tion events originate from electrons within the second con- 
duction band. The number of events originating in the first 
conducticm band is found to be higher for the numerical 
model than for the model using the Quade formula at the 
field of 250 kV/cm, while for higher fields this difference 
vanishes. His slightly different behavior of the two models 
at low fields is caused by the fact that the threshold energy 
for impact ioniation is about 132 eV in the Quade model 
for the particle originating in the first conduction band and 
2.11 eV for the second conduction band, while for the nu- 
merical model the threshold energy value is not well defined. 



CrMrgy(»V) 

FIG. 6. The electron number density function in silioon, defined as the 
product of the density-of^ststes function end the electron distribution func- 
tioQ, It an applied electric field of 500 kV/cm, calculiied using the 
h-depeadent ouinerical model and the Quade model. 


Instead, in the numerical model ionization is possible at all 
states with energy greater than the band gap which are found 
to satisfy the momentum and energy conservation conditions. 
Therefore, there is t» sharp threshold energy, other than the 
band gap, below which ionization cannot occur and above 
which ionization does occur. Thus, the contribution of the 
first conduction band at low fields can be expected to be 
higher when the numerical ntodel is used. 

The effect of the second conduction band on the impact- 
ionization rate in gallium arsenide is even more important 
than in silicon. The impact-ionization transition rate for elec- 
trons in the first conduction band is much lower than the rate 
within the second conduction band in both models as can be 
seen frcmi Figs. 2(a) and 2(b) and from Ref. 16. As a result, 
'virtually all the electrons which impact ionize originate 
within the second conduction band. The effect is even stron- 
ger using the Quade formulation than the numerical formu- 
lation. Though the value of F used in the analytical formu- 
lation is unknown a priori, i.e., without first comparing the 
calculations to experiment, it should be noted that F can 
never be greater than 1. Therefore, given that the transition 
rate of electrons within the first conduction band is several 
orders of magnitude lower than from the second conduction 
band with 1, as shown by Figs. 2(a) and 2(b), it is clear 
that the first conduction band plays no significant role in 
impact ionization in bulk GaAs for reasonable values of F. 
Hence, in either the numerical or analytical k-dependent 
modek the predominant source of ionizing electrons in bulk 
GaAs is carriers within the second conduction band. 

The election number density function defined as the 
density-of-states function multiplied by the electron distribu- 
tion fimetion, is shown in Fig. 6 iot silicon at an applied 
electric field of 500 kV/cm for both the numerical 
k-dqiendent model and the Quade model. In both cases, the 
number of electrons which survive to high energies is rela- 
tively high. The second peak in each of the figures corre- 
sponds to the electron population in the second conduction 
band and is slightly larger for the Quade model. Both models 
essentially confirm previous results^^'^’ that the threshold in 
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FIG. 7. The electron number density function in fsllium arsenide, defined as 
the product of the density-of-ststes function and the electron distribution 
function, at an applied electric field of 500 kV/cm, calculated using the 
k'dependent numerical model and the Quade model. 

silicon is exceedingly soft and that the electrons can reach 
energies substantially higher than the band gap before they 
impact ionize. 

The electron number density functions in gallium, ars- 
enide are plotted in Fig. 7 calculated using the k-dependent 
numerical and Quade models for the applied electric field of 
500 kV/cm, In both cases the distribution functions are simi-/ 
lar and the effect of the second conduction band (the second 
peak at higher energies) is again slightly larger in the Quade 
model. Notice that the number distribution remains substan- 
tially large at high carrier energies indicating an exceedingly 
soft effective threshold for impact ionization in GaAs. Again 
these results are consistent with previously reported 
observations.^ 

The energy-dependent impact-ionization transition rate, 
derived from the k-dependent rate in a manner described in 
Ref. 15, is plotted in Fig. 8 for the numerical k-dependent 
formulation and the Quade formula, respectively, for silicon. 



FIG. 8. Calculated impact^kmizatton transition rate in tfliooa as a function 
of electron energy measured from the conduction-band minimum. The 
energy-dependent rate is calculited by performing a further integration of 
the k-dependent rate following the approach explained in Ref. 15 for the 
numerical and Quade models; the results of Kane (Ref. 40) and CSrtier tt 0 L 
(Ref. 41) and Kamakura et al (Ref. 19) are also shown. 


FIG. 9. Cklculsted impact-kmizslion transition rate in gallium arsenide as a 
function of election energy measured from the conduction-band minimum 
■sing the k-dependent numerical model and the Quade model. The energy- 
dependent rate is calculated by performing t further iniegralion of the 
k-dependent rate following the a pproach explained in Ref. 15. The energy- 
dependent rate obuined by Stobbe and oo-workers (Ref. 21) is also in- 
cluded. 

The results of Kane,^ Kamakura et and the recently 
published ionization transition rate of Cartier et al.*' deter- 
mined from a comparison of Monte Carlo simulations and 
soft-x-ray photoemission spectroscopy (XPS) measurements 
are shown as well. The numerically calculated transition rate 
is almost identical with that of Cartier et al. up to an energy 
of 4 eV, and teems to be closer to the results of Kane in the 
high-energy region. Recently, also, Kunikiyo era/.‘* have 
performed similar calculations of the sUioon interband 
impact-ionization tninsition rate. Their calculations are in ex- 
cellent agreement as well with the transition rate of Cartier 
et al.** The transition rate obtained on the basis of the Quade 
formula agrees well with the numerical model between 2 and 
4 eV, while a lower rate is predicted in the low- and high- 
energy regions. 

The energy-dependent impact-ionization transition rate 
for gallium arsenide, calculated using the numerical ap- 
proach of Wang and Brennan'* and the Quade formula, re- 
spectively, arc shown in Fig. 9. Stobbe snd co-workers * 
have also recently calculated the impact-ionization transition 
rate in bulk OaAs. Their calculations, based on the fitting 
formula quoted in Ref. 21, are included in Fig. 9 for com- 
parison. Inspection of Fig. 9 shows that the maxima of the 
transition rate, based on both the numerical model of Wang 
and Brennan*® and the Quade model, occurs near 3 eV in 
OaAs. The maximum region is more pronounced for the 
Quade formulation than in the case of Wang and Brennan’s 
numerical model. The origin of the maximum can be possi- 
bly understood within the analytical formulation, since the 
transition rate involving each valley can be individually 
evaluated. The maximum in the transition rate is due to the 
collision of an electron within the second conduction band 
with an electron in the heavy-hole band, after which the elec- 
trons occupy stales either both in the F valley (the lower 
peak below 3.0 eV) or one is In the F valley and the other in 
the L valley (the higher peak close to 3.4 eV). In both cases, 
a transition of the initiating electron from the second conduc- 
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tton band to the first one with only a slight change of the k 
vector (a near-vertical transition, as explained in Ref. 42) 
contributes to the final values of the transition rate. Near- 
vertical transitions have a significantly greater rate owing to 
the fact that the denominator in Eq. (3) becomes very small 
under these conditions.^^ This region may be somewhat un- 
realistically pronounced in the Quade formulation due to the 
fact that the value of the overlap integrals is taken as an 
averaged constant. The actual values of the overlap integrals 
in this region may be quite small thereby reducing the rate if 
calculated. It is much more difficult to isolate the physical 
origin of the maximum in the numerical model. In this case, 
a similar classification would be questionable, mainly due to 
the uncertainty in how the *‘valleys” are defined in the real 
band structure. Inspection of Fig. 9 also shows that the study 
of Stobbe and oo-workers^^ predicts a higher ionization tran- 
sition rate from the results presented here. The cause of this 
discrepancy may be due to the inclusion of a dependent 
dielectric function in Stobbe and co-workers’ formulation as 
well as differences in the energy band structure utilized in 
the two calculations. In Stobbe and co-workers* work the 
band structure is calculated using an empirical pseudopoten- 
tial method in contrast to the k<p method adopted in the 
present work. 

To further investigate the impact-ionization transition 
rate formulations, the energy dependence of the quantum 
yield was also calculated for both the numerical and Quade 
models. The quantum yield is defined as the average number 
of impact-ionization events caused by a high-energy injected 
electron until its kinetic energy relaxes below the ionization 
threshold through scattering and/or ionization events. While 
the impact-ionization coefficient strongly depends on the en- 
ergy distribution function and therefore mainly on the inelas- 
tic electron-phonon scattering rate, the quantum yield is al- 
most a linear function of the ratio of the ionization transition 
rale to the electron-phonon scattering rate.^' The quantum 
yield provides a useful means of assessing the relative hard- 
ness or softness of the impact-ionization threshold and en- 
ables an additional, independent evaluation of the appropri- 
ateness of the model used. 

To the authors* knowledge, the only presently existing 
experimentally based data for the quantum yield apply to 
bulk silicon. Calculations of the quantum yield in silicon are 
presented in Fig. 10 using both the numerical and Quade 
models. The data of Cartier et derived from a compari- 
son of Monte Carlo simulations and XPS measurements, 
which were shown to be in good agreement with other 
measurements,^ ‘ are also plotted in Fig. 10 for comparison. 
While fairly good agreement over the full energy range be- 
tween the model based on the numerically calculated ioniza- 
tion rate and the data from Ref. 41 can be seen, the Quade 
formula predicts a quantum yield which is substantially 
lower in the low-energy region. Nevertheless, the Quade for- 
mula results correspond fairly well to the other data above 
2.5 eV. These results can be understood on the basis of the 
threshold energy for impact ionization. In the case of the 
Quade formula the threshold energy is well defined and is 
substantially higher than the value of the gap. Subsequently 
the ionization-transition rate in this region is lower for the 



FIO. 10. GdcuUted ipiuitum yield as • fuactioa of electron energy in till' 
coo metsured from the oonduction-hend minimum, ohuined using Uie 
k-depeadenl anmerictl model nod Ute Quade model. Results of Cartier et at 
(Ref. 41) arc abo shown. 

Quade model, as can be aeen in Fig. 8, which results in a 
lower quantum yield. Conversely, in the numerical model, 
the effective threshold is simply the energy gap and subse- 
quently ionization events can possibly occur at low energies 
in this model leading to a higher quantum yield. 

Although no experimental data for bulk GaAs are pres- 
ently available, it is nevertheless instructive to compare the 
calculated quantum yield dependencies on energy using the 
numerical and Quade models. The calculated quantum yield 
versus electron energy in GaAs using the numerical and 
Quade models is shown in Fig. 11. Comparison of the calcu- 
lations reveals good correspondence between these two mod- 
els in the energy range from 2.S to eV, while at higher 
energies the Quade model [with the particular choice of the 
prefactor F appearing in Eq. (7), as described above] predicts 
a lower quantum yield than the numerical model. At lower 
energies the Quade model predicts a higher quantum yield 
than the numerical model. This discrepancy can be under- 
stood as follows. An identical phonon scattering rate is used 
in both models of the ionization rate. As stated by Cartier 
et the quantum yield depends linearly upon the ratio of 
the impact-ionization transition rate to the phonon scattering 



FIO. 11. Ctlculaied quantum yield as a function of electron energy in gal- 
lium aneoide measured from the conduction-band minimum, obuined using 
the k-dependent numerical model and the Quade model. 
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rate. Subsequently, the quantum yield within each model, 
given that the scattering rales are identical, depends then 
simply upon the ratio of the transition rates and will follow 

the same trends as shown in Fig. 9. 

The value of the overlap integrals used in the evaluation 
of the transition rate given by Eqs. (2) and (3) is evaluated 
directly in the numerical approach. As mentioned above, this 
requires extensive numerical computation. By approximating 
these overlap integrals as a constant, the factor F in Eq. (7), 
a new single parameter is introduced greatly redudng the 
numerical complexity required to evaluate the transition rate 
but at the expense of a new unknown parameter. In the cal- 
culations presented above, F is adjusted such that the ana- 
lytical model yields satisfactory agreement with the experi- 
mental data. To check to see if the determined values of F 
are reasonable, it is useful to compare these values to those 
determined using the sum rule quoted by Ridley. Using the 
sum rule, the squared overlap integral between the first con- 
duction band and the heavy-hole band in gallium arsenide is 
determined to be 0.38. More elaborate calculations made by 
Burt et al.*^’** using wave functions obtained from a k*p 
calculation show that this value is which is in fair 

agreement with the results presented here for GaAs, 
F=0.005 for the second conduction band. In the case of 
silicon, the value F= 0.24 obtained for the impact ionization 
originating in both the first and second conduction bands, is 
higher than the value of F-0.05,"* or F=0.01. quoted in 
Refs. 46 and 47, obtained by theoretical evaluations. It 
should be noted, though, that a discrepancy between the ex- 
perimentally determined recombination coefficients for Au- 
ger recombination, and the theoretically obtained ones, 
which also contain the same overlap integrals, was observed 
in bulk silicon. This disaepancy is thought to be possibly 
due to an underestimation of the overlap integrals by roughly 
an order of magnitude.*® 


IV. CONCLUSIONS 

in this article the electron interband impact-ioiiization 
rate calculated by an ensemble Monte Carlo simulation has 
been presented for both bulk silicon and gallium arsenide 
with the expressed purpose of comparing different formula- 
tions of the transition rate. Specifically, three different Rai- 
ments of the transition rate are examined: the traditional 
Keldysh formula, a new k-dependent BMlytical formulation 
first derived by Quade and co-workers, and a more exact, 
numerical method of Wang and Brennan. ’ The numerical 
method of Wang and Brennan'* '® conUins no readily adjust- 
able parameters and is similar to other numerical tRhniquM 
that arc currently attaining prominence. The numerical tech* 
nique thus serves as a standard to which the Quade result has 
been compared. Unfortunately, the completely numental 
technique, although fundamentally based, requires extensive 
numerical computation making it somewhat unattractive. ^ 
is well known, the Keldysh formula, although commonly 
used in Monte Carlo calculations of the impact-ionization 
rale owing to its relative simplicity and computational effi- 
ciency, suffers from several limitations as well. Chief among 
these is the fact that it is derived assuming parabolic energy 


bands for all of the participating carriers, and includes no k 
dependence in the rate. Additionally, due to its reliance on 
two adjustable parameters, more than one set of parameters 
can be used within the Keldysh formula to produce satisfac- 
tory agreement with experiment. As a result, the Keldysh 
formula cannot by itself predict whether a material has a soft 
or hard ionization threshold. Therefore, there is a need to 
develop a more computationally efficient transition rate for- 
mulation than the direct, numerical technique, yet retain the 
k dfpfnd>n«. and deUib of the band structure in its formu- 
lation. In thb way, the physical future of the ionization pro- 
cess can be preserved (soft or hard threshold, etc.) while 
maintaining compuUtional efficiency. The Quade formula, 
although it contains one adjusUble parameter corresponding 
physically to the squared overlap integrals of the periodic 
parts of the Bloch functions, includes the full details of the 
band structure for the initiating carrier. In this article, we 
have provided the first comparison of thb formula to the 
numerical tedinique to determine if this aprfoach yields 
qualitatively similar results while providing far better com- 
putational efficiency. 

The calculations presented here show that the model 
based on the Quade formulation produces similar results to 
the completely numerical calculations for some quantities, 
primarily the ionization rate as a function of field. In con- 
trast, to the Keldysh formula, where no k-vector dependence 
of the impact-ionization transition rate is included, the tran- 
sition rate in the k-dependent numerical and analytical 
(Quade formula) cases b shown to be strongly dependent on 
the initiating electron k vector. It b also shown that both 
k-dependent modeb predict that the second conduction band 
dominates the ionization process in GaAs and silicon; how- 
ever, the quantum yield predicted using the analytical model 
of Quade and co-workers“ does not agree cItRly with the 
available experimentally derived data of Cartier et al. or 
that predicted using the numerical model of Wang and 
Brennan.'* '® As such the Quade formula, although far more 
complete than the Keldysh formula, may not be a completely 
reliable rqilacement for the numerical models in many situ- 
ations. 

Further comparison of the Quade and numerical models 
reveab that they predict that the effective threshold for im- 
pact ionization b rehtively soft, implying that the Rfriers 
drift to substantially high energies before suffering an ioniza- 
tion event. It b shown that the average energy of the elec- 
trons b high when compared to that expected for a hard 
Ihicdiold model. Tlius, many of the Important qualitative re- 
sulb obttined by adopting the computationally intensive nu- 
merical approach, Le., the dominance of the second conduc- 
tion band, the nature of the ionization threshold, and the 
importance of the k dependence, can be recovered if the 
analytical k-dependent fwmulation of Quade and 
oo-workeis^ b used. Use of the Quade formula may not 
always be fully justifiable, however, since it contains some 
ad hoc parametrization which cannot be ascertained from 
first principles and it apparently does not reproduce the quan- 
tum yield data. Additionally, it b possible that the Q^de 
formula may overestimate the impact-ionization transition 
rate at those points in the k space where near-vertical transi- 
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tions are allowed due to the failure to include the actual 
values of the overlap integrals for these transitions. 
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